Report

Models for small area data with applications in health care Nicky Best Department of Epidemiology and Biostatistics School of Public Health, Imperial College London http://www.bias-project.org.uk http://www1.imperial.ac.uk/medicine/people/n.best/ Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data A brief history of disease mapping Health indicator maps have a long history in epidemiology and public health • Spot maps: – Yellow fever pandemic New York (Seaman, 1798) – Cholera and the Broad Street Pump (Snow, 1854) Spot map of cholera cases (Snow, 1854) A brief history of disease mapping Health indicator maps have a long history in epidemiology and public health • Spot maps: – Yellow fever pandemic New York (Seaman, 1798) – Cholera and the Broad Street Pump (Snow, 1854) • Chloropeth maps: – Geographical distribution of mortality from heart disease, cancer and TB in England & Wales (Haviland, 1878) – Cancer mortality by county in England & Wales, adjusted for age and sex (Stocks, 1936, 1937, 1939) Female cancer 1851-60 (Haviland 1878) Female lung cancer SMR 1921-30 (Stocks, 1939) A brief history of disease mapping Health indicator maps have a long history in epidemiology and public health • Spot maps: – Yellow fever pandemic New York (Seaman, 1798) – Cholera and the Broad Street Pump (Snow, 1854) • Chloropeth maps: – Geographical distribution of heart disease, cancer and TB in England & Wales (Haviland, 1878) – Cancer rates by county in England & Wales, adjusted for age and sex (Stocks, 1936, 1937, 1939) • National and international disease atlases, e.g – Atlas of Cancer Incidence in England & Wales 1968-85 (Swerdlow & dos Santos Silva, 1993) – Atlas of Mortality in Europe 1980/81 & 1990/91 (WHO, 1997) Female lung cancer incidence 1968-85 (Swerdlow and dos Santos Silva, 1993) Age-standardised mortality from IHD, 1980-81 (WHO) Recent developments in disease mapping • Development of Geographical Information Systems (GIS) – Geographically indexed relational database – Computer program to map and analyze spatial data • Increasing availability of geo-referenced data – Ability to geocode, use GPS – Disease outcomes, demographics, environmental quality, health services • Development of statistical methods – Sophisticated techniques for separating signal from noise – Ability to account for spatial (and temporal) dependence – Methods for cluster detection and classification of areas Interest in mapping health events at small-area scale Small area health data in the UK • Administrative geography in UK includes – Postcodes (10-15 households) – Census Output Areas (COA; ~300 people) – Electoral wards (~500 to 2000 people) – Local authority districts, Health authority districts (10’s of thousands) • Postcoded data on mortality, births/still births, congenital anomalies, cancer incidence, hospital admissions • Population and socio-economic indicators from Census (COA) • Increasing availability of modelled environmental data at fine geographical resolution (grids) • Limited access to geographical identifiers for certain individual-level cohorts (e.g. Millenium Cohort, British Household Panel Survey) and health surveys (e.g. Health Survey for England) Small area health data in Spain • Administrative geography in Spain is divided into: – 17 regions – 52 provinces – ~8000 municipalities, ranging from small villages to large cities – Census tracts (finer sub-division in large cities) • Geocoded (place of residence) data on births, mortality (national), cancer incidence (regional; ~26% population), hospital discharge administrative data (national; public hospitals) • Small area (municipality) data on population and socioeconomic indicators from Census Examples of recent disease atlases and health-related maps for Spain Atlas of cancer mortality and other causes of death in Spain 1978-1992 (López-Abente et al., 1996) Maps Ageadjusted Rates and Standardised Mortality Ratios (SMR) at province level Atlas of cancer mortality at municipality level in Spain 1989-1998 (López-Abente et al., 2007) Maps (Bayesian) smoothed relative risks of mortality and probability of excess risk, at municipality level Also produced maps of mortality from selected causes other than cancer, e.g. Influenza …… + contextual maps of socioeconomic variables and environmental hazards Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data Why map small area disease rates? • Interest in mapping geographical variations in health outcomes a the small area scale – Highlight sources of heterogeneity and spatial patterns – Suggest public health determinants or aetiological clues • Small scale – UK: electoral ward or census output area – SPAIN: municipality – less susceptible to ecological (aggregation) bias – better able to detect highly localised effects Why smooth small area disease rates? • Typically dealing with rare events in small areas Ai Yi is the observed count of disease in area Ai Ei is the expected count based on population size, adjusted for age, sex, other strata …., • Relative risk usually estimated by SMRi = Yi / Ei • Standard practice is to map SMRs BUT sparse data need more sophisticated statistical analysis techniques Why smooth small area disease rates? • SMR represents estimate of ‘true’ (underlying) risk in an area, Ri, i.e. Ri = SMRi • Statistical uncertainty about estimate based on assuming Poisson sampling variation for data Yi ~ Poisson(Ri Ei) • SE(Ri) = SE(SMRi) 1 / Ei – SMRi very imprecise for rare diseases and small populations – precision can vary widely between areas Why smooth small area disease rates? • SMRi in each area is estimated independently – ignores possible spatial correlation between disease risk in nearby areas due to possible dependence on spatially varying risk factors – leads to problems of multiple significance testing Map of SMR of adult leukaemia in West Midlands Region, England 1974-86 (Olsen, Martuzzi and Elliott, BMJ 1996;313:863-866). Is the variability real or simply reflecting unequal expected counts ? Have the red highlighted areas truly got a raised relative risk? Methods for smoothing disease maps • These problems may be addressed by spatial smoothing of the raw data • Idea is to “borrow information” for neighbouring areas to produce better (more stable, less noisy) estimates of the risk in each area • Similar principle to scatter plot smoothers, moving average smoothers…. • Many methods available Methods for smoothing disease maps • Ad hoc, local smoothing algorithms e.g. spatial moving averages, headbanging algorithm + quick and simple to implement – can be very sensitive to ad hoc choice of weights etc. – no uncertainty estimates (standard errors) • Trend surface analysis e.g. kriging, polynomial/spline smoothing + estimation of ‘smoothing parameters’ based on tradeoff between fit and smoothness – can be sensitive to choice of penalty for trade-off + standard errors usually available Methods for smoothing disease maps • Random effects models e.g. empirical Bayes, hierarchical Bayes + data-based estimation of model parameters that control smoothing + full power of statistical modelling available: standard errors, prediction, probability calculations, inclusion of covariates – more complex to understand and implement Bayesian Approach • Use probability model to obtain smoothed risk estimate Ri in area i that is a compromise (weighted average) of – observed area-level risk ratio (Yi/ Ei) – local or regional mean relative risk (m) • Weights depend on the precision of the SMR ( 1 / Ei) in area i and the variability (heterogeneity) of the true risks across areas (v) local or regional mean relative risk (m) • Aim is to estimate posterior probability distribution of the unknown model parameters (Ri, m, v) conditional on the data (Yi/ Ei) Bayesian disease mapping model • Typical Bayesian disease mapping model: Yi ~ Poisson(Ri Ei), log (Ri) ~ Normal (m, v) • Hierarchical Bayesian model also requires specification of a (prior) probability distribution for m and v – These are often taken to be ‘non-informative’ • Empirical Bayes involves 2-step process: – Estimate m and v empirically from observed data – Ignore uncertainty in estimates of m and v and plug these values into the Bayesian model above Software • Estimation of Bayesian hierarchical models requires computationally intensive simulation methods (MCMC) – Implemented in free WinBUGS and GeoBUGS software: www.mrc-bsu.cam.ac.uk/bugs • Free software INLA (Rue et al, 2008) implements fast approximation: www.r-inla.org • Empirical Bayes smoothing implemented in Rapid Inquiry Facility (RIF): www.sahsu.org/sahsu_studies.php#RIF Map of occurrences of adult leukaemia in West Midlands Region, England 1974-86 (A) unsmoothed SMR (B) smoothed by Bayesian methods (Olsen, Martuzzi and Elliott, BMJ 1996;313:863-866) 0 5 Expected count Expected 0count 2 4 Comparison of estimation methods relative risk 0 2 4 6 8 SMR Hierarchical Bayes Empirical Bayes 5 10 area 15 20 Including spatial dependence in disease risk • Ri are typically spatially correlated because they reflect, in part, spatially varying risk factors Incorporation of spatial dependence in the distribution of the Ri’s Conditional Autoregressive (CAR) model log (Ri ) ~ Normal (mi , vi) mi = k Rk / ni = average risk in neighbouring areas vi = v / ni → variance inversely proportional to number of neighbours Besag, York, Mollie (1991) Annals of the Institute of Statistics and Mathematics, 43: 1-59 Childhood leukaemia incidence in London, 1986-1998 (summary)means for RR Non-spatial smoothing (posterior mean Ri) N SMR Raw data (SMR) (418) < N (7) 0.5 0.5 - RR (30) 0.7 0.5 <0.5 (55) 0.9 - 0.7 0.9 (0) < (39) 0.5 0.5-0.7 (71) 1.1 0.7 0.7-0.9 (109) 1.4 - (192) 1.1 0.7 1.4 RR 0.9 2.0 N 0.9-1.1 (264) (183) 0.9 >= - 2.01.1 1.1 1.1-1.4 (251) 1.4 1.4-2.0 (96) >2.0 (31) >= 20.0km 1.4 2.0 2.0 20.0km Spatial smoothing (posterior mean Ri) 20.0km 20.0km Mapping uncertainty • Mapping the mean posterior value of Ri does not make full use of the posterior distribution • Map posterior SD • Map Probability (Ri > 1) Note – this is not the same as a classical p-value 0.5 0.5 1.0 1.0 1.5 1.5 Relative Risk 2.0 2.0 Relative Risk, Ri 2.5 2.5 Posterior SD of relative risk estimates Posterior mean relative risk Posterior sd of relative risk RR values for sd(0) < N (39) 0.5 0.5 - 0.7 (0) < 0 (324) 0 0 N (192) 0.7 - 0.9 (418) (264) 0.9 - 1.1 (91) 0. (251) 1.1 - 1.4 (20) 0. (5) 1.0 (96) (31) >= 1.4 - 2.0 (15) >= 2.0 RR SD (0) < 0.5 values for sd <0.5 (39) N 0.5 0.5-0.7 0.2 <0.2 0.7 (324) 0.2 0.4 0.2-0.4 0.7 0.7-0.9 0.9 0.9 0.9-1.1 1.1 (91) 0.6 - 0.8 0.6-0.8 1.1 1.1-1.4 1.4 (20) 0.8 - 1.0 0.8-1.0 (5) 1.0 1.2 1.0-1.2 (192) 20.0km (0) < (264) (251) 1.4 1.4-2.0 (96) >2.0 (31) >= 2.0 2.0 20.0km (418) 0.4 0.6 0.4-0.6 >1.21.2 (15) >= Posterior probability that relative risk > 1 Posterior mean relative risk Posterior probability that relative risk > 1 values for pp(0) < RR N (39) N (163) 0.5 0.5 - (321) 0.7 (192) 0.7 - 0.9 (250) (264) 0.9 - 1.1 (139) (251) 1.1 - 1.4 (96) (31) >= 1.4 - 2.0 2.0 RR Prob (0) < <0.50.5 values for pp (39) 0.5 0.5-0.7 20.0kmN (163) < 0.25 <0.25 0.7 (321) 0.25 - 0.9 0.9 0.9-1.1 1.1 (250) 1.1 1.1-1.4 1.4 (139) >= 0.75 (264) (251) 1.4 1.4-2.0 (96) >2.0 (31) >= 2.0 20.0km 2.0 0.5 0.25-0.50 0.7 0.7-0.9 (192) 0.5 - 0.75 0.50-0.75 >0.75 Atlas of cancer mortality at municipality level in Spain 1989-1998 (López-Abente et al., 2007) Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data Classifying areas with excess risk • Richardson et al (2004): Simulation study investigating use of posterior probabilities in disease mapping studies • Classify an area as having an elevated risk if [Prob (Ri > 1)] > 0.8 • High specificity (false detection < 10%) • Sensitivity 60%-95% for Ei of 5-20 and true Ri of 1.5-3.0 Childhood leukaemia in London Posterior prob(RR>1) Posterior mean RR probability of(0) RR < greater 0.5 than 1.0 RR N N (39) 0.5 - 0.7 (192) 0.7 - 0.9 (264) 0.9 - 1.1 (251) 1.1 - 1.4 (96) (31) >= 1.4 - 2.0 2.0 (0) < <0.50.5 (39) 0.5 0.5-0.7 0.7 0.7-0.9 (192) 0.7 0.9 probability of RR greater than 1.0 0.9 - 0.9-1.1 1.1 1.1 1.1-1.4 1.4 (264) 20.0km N (251) 1.4 1.4-2.0 (96) >2.0 (31) >= 2.0 2.0 20.0km (772) <0.8< 0.8 (101) >= >0.8 0.8 Comparison of SaTScan and Bayesian classification rule SaTScan (Kuldorff, www.satscan.org): Location of most likely cluster Bayesian: probability of RR greater than 1.0 of excess risk Probability N probability of RR greater than 1.0 Most likely cluster; p<0.001 2nd most likelyN cluster; p = 0.2 20.0km (772) < 0.8 < 0.8 (101) >= 0.8 ≥ 0.8 Summarising geographic variation • Often interested in providing overall summary measure of variability between areas, e.g. – to compare variability of different outcomes – to quantify how much variation can be explained by covariates • Percentile Ratio: Ratio of outcomes (relative risks) in areas ranked at the qth and (100-q) th percentiles – e.g. 90th Percentile Ratio, PR90 = R95%/R5% – Posterior distribution of PR90 easily calculated from MCMC output Relative survival from colon cancer, England Data • Survival/censoring times for all 7007 cases of colon cancer diagnosed in England in 1995 and followed for 5 years (provided by B Rachet, LSHTM) • Covariates: sex, age at diagnosis, clinical stage, deprivation score, Health Authority (95 area, 1-300 cases per HA) • Population mortality rates by age and sex for England and Wales, 19851995. Questions of interest • Is there evidence of differences between Health Authorities in relative survival that may indicate differences in effectiveness of care received? – Relative survival measures difference between age/sex-adjusted mortality rate in general population and in patients with disease of interest • How do these geographical differences change when we adjust for socioeconomic deprivation and clinical stage of cancer? Relative survival from colon cancer, England ykit ~ Poisson(mkit) (subject k, area i, time interval t) log(mkit – Ekit) = log nkit + at + bxki + Hi Standard model for relative survival Area spatial effect Relative survival from colon cancer, England ykit ~ Poisson(mkit) (subject k, area i, time interval t) log(mkit – Ekit) = log nkit + at + bxki + Hi Area spatial effect Standard model for relative survival adjustment relative Without exces s hazard (model 6) for deprivation and clinical stage (7)relative < 0.85 (31) After adjustment exces s hazard (model 3)for deprivation and clinical stage 0.95 0.85 - (7) (16 N N Relative (28) 0.95 excess (9) 1.05 mortality (36 1.05 (25 1.15 (11 (20) >= 1.15 relative exces s hazard (model 6) (7) < 0.85 <0.85 (31) 0.85 0.85-0.95 0.95 (28) 0.95 0.95-1.05 1.05 N (9) 1.05 - 1.15 1.05-1.15 PR90 = 1.95 (95 % CI 1.62-2.38) (20) >= 1.15 ≥1.15 PR90 = 1.83 (95 % CI 1.54-2.24) Ranking and classifying extreme areas • Interest in ranking areas for e.g. policy evaluation, ‘performance’ monitoring • Rank of a point estimate is highly unreliable – Would like to measure uncertainty about rank • Straightforward to calculate posterior distribution of ranks (or any function of parameters) using MCMC – Obtain interval estimates for ranks • Can also calculate posterior probability that each area is ranked above a particular percentile Rank (posterior mean and 95% CI) of the 95 Health Authorities Without adjustment for deprivation and clinical stage caterpillar plot: rank.HA After adjustment for deprivation and clinical stage caterpillar plot: rank.HA Upper quartile 0.0 25.0 50.0 Rank 75.0 Upper quartile 0.0 25.0 50.0 Rank 75.0 Posterior probability that HA is ranked in top 5% After adjustment for deprivation and clinical stage Without adjustment for deprivation and clinical stage (samples)means for p.bottom5 (samples)means (20) < 0.0 for p.bottom5 N N (64) (6) (samples)means for p.bottom5 0.0 0.1 - 0.0 (4) 0.2 0.0-0.1 (1) >= 0.5 0.1-0.2 0.2-0.5 >0.5 (20) < 0.2 0.0 (64) N 0.1 0.0 - (6) 0.1 - 0.2 (4) 0.2 - 0.5 (1) >= 0.5 0.1 0.5 200.0km 200.0km 200.0km Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data Joint spatial variation in risk of multiple diseases Disease 1 Disease 2 RR SMR 0.7 0.7 SMR 0.7 Specific component 1 RR 0.7 1.0 1.5 1.0 1.5 SMR 0.7 Shared component 1.0 1.0 1.0 1.5 Specific component 2 RR Knorr-Held and1.5Best (2001) RR 0.7 1.0 0.7 1.0 1.5 1.5 1.5 Statistical model Y1i ~ Poisson(R1i E1i); log R1i = Si + U1i Y2i ~ Poisson(R2i E2i); log R2i = Si + U2i Si ~ spatial model (shared component of risk) U1i ~ spatial model (component of risk specific to disease 1) U2i ~ spatial model (component of risk specific to disease 2) • Extends to >2 diseases (Tzala and Best, 2006) • Extends to shared variations in space and time (Richardson et al, 2006) © Imperial College London Joint variation in COPD and lung cancer in GB values for SMR1 COPD SMR N values for SMR2 (42 ) < 0.66 (10 0) 0.66 - 0 .8 (13 9) 0.8 - 1 .0 (12 1) 1.0 - (40 ) 1 .25 - (17 ) >= Lung cancer SMR N (11 ) < 0.6 (10 1) 0.6 (20 1) 0.8 1.25 (10 3) 1.0 1.5 (33 ) (10 ) >= 1.5 200.0km 200.0km Best and Hansell (2009) 1 .25 1 Modelled risk estimates values for shared Shared risk N Shared risk interpreted as mainly reflecting geographical variations in community-level smoking behaviour (16 8) < values for log.COPD.spec.MALE.model2 0.9 (40 ) 0 .9 - (32 ) 0 .95 - (34 ) 1 .0 - (25 ) 1 .05 - (28 ) 1 .1 - (13 2) >= COPD specific risk 0 .9 5 1.0 N 1 .1 5 1.15 -0.1 - (111) COPD specific risk interpreted as reflecting smoking-adjusted variations in COPD mortality 200.0km 200.0km (39) (71) -0.05 - 1 .0 5 1.1 (80) < -0.1 0.0 - (89) 0.05 - (36) 0.1 - (33) >= 0.15 Joint variation in relative survival of colon and breast cancer by English Health Authority • Shared spatial patterns of relative survival may reflect variations in effectiveness of health care system Observed 5-year relative survival: Breast MLE estimate of relative survival: breast Observed 5-year relative estimate of relative survival: colon survival: Colon 65.0 - 70.0 (3) < MLE 65.0 (13) N N (36) 70.0 - 75.0 (3) 65.0 < <65% (38) 75.0 - 80.0 (13) 70.0 65%65.0 to -70% (5) >= 80.0 MLE estimate of relative survival: colon (2) < 20.0 < 20% (10) - 30.0 20%20.0 to 30% N (36) 75.0 70%70.0 to -75% (31) - 40.0 30%30.0 to 40% (38) 80.0 75%75.0 to -80% (36) - 50.0 40%40.0 to 50% (5) >= 80.0 >80% (16) >= 50.0 >50% 200.0km 200.0km Difference in relative survival in each HA compared to England as a whole Posterior Prob that shared difference > 0 Shared difference shared differential in relative s urvival (3) < -3.0 of shared > 0.0 probability N N (13) -3.0 - -1.5 (17) < -1.5 0.2 <(63) 0.2 1.5 (13) -3.0 -1.5 -15% to --30% (15) (64) 0.2 3.0 0.8 (63) -1.5 1.5 -15% to -15% (14) >= 0.8 > 0.8 < -3.0 -30% prob<(3) that relative excess hazard N 1.5-–0.2 0.8 (1) >= (11) < (73) (11) > 3.0 (15) to 1.530% 3.0 15% (1) >= 3.0 >30% 200.0km 200.0km Difference in relative survival in each HA compared to England as a whole Difference specific to colon cancer Difference specific to breast cancer breast specific differential in rel. surv (0) < -3.0 colon s pecific differential in rel surv. (1) N -3.0 - -1.5 N < -3.0 differential in relative s urvival <(0)shared -30% (94) (3)-30% <- -3.0 < -1.5 1.5 (1) -3.0 -1.5 -15% to- -30% (0) (13) 1.5 - (94) -1.5 1.5 -15% to 15% (0) >= (63) 3.0 N -3.0 - -1.5 -15%3.0 to -30% 1.5 -15%-1.5 to 15% (0) 1.5 3.0 15% to -30% (15) to 1.530% 3.0 15% (0) >= >30% (1) >= >30% 3.0 200.0km 3.0 200.0km Diet-related cancers in Greece Cancer-specific spatial residuals Spatial common factor 1 2 3 4 5 values for p oesophagus (13) <-0.1184 values for p N N(10) -0.1184 - -0. 02204 (6) -0.11 (8) -0.02204 - 0.03013 (1) -0.01 (8) 0.03013 - 0.1332 (7) 0.030 (12) >= 0.1332 (19) >= (10) -0.005441 - -0.002873 N stomach values for p (11) <-0.005441 (10) -0.002873 - -1.145E-4 (18) <-0 (10) -1.145E-4 - 0.001841 (10) >=0.001841 200.0km values for p N colorectal 200.0km values for p (14) <-0.1234 N(10) -0.1234 - -0. 01909 pancreas (6) <-0.1 (4) -0.12 (7) -0.01909 - 0. 0284 (2) -0.01 (17) 0.0284 - 0. 1359 (15) 0.0 (3) >= 0.1359 (24) >= 200.0km 200.0km 200.0km values for p N (0) <-0.1234 values for p prostate (18) -0.1234 - -0. 01909 N Cut-points based on quintiles of distribution of factor values and of residuals across all cancers (10) <-0.11 bladder (9) -0.1184 (30) -0.01909 - 0.0284 (14) -0.022 (3) 0.0284 - 0.1359 (15) 0.0301 (0) >= 0.1359 (3) >= 0.1 Tzala and Best (2006) 200.0km 200.0km Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data Extensions of disease mapping to space time modelling Noisy data in each area Noise model: Poisson/Binomial Latent structure: Space + Time + (Residuals) + joint Bayesian estimation Inference Basic space-time model set-up Yit ~ Poisson(Rit Eit); log Rit = Si + Tt + Uit Si ~ spatial CAR model (common spatial pattern) Tt ~ random walk (RW) model (common temporal trend) Uit ~ Normal(0, v) (space-time residual reflecting idiosyncratic variation) • Extends to shared variations of 2 outcomes in space and time © Imperial College London Space-time variations in Male and Female lung cancer incidence (Richardson et al, 2006) • Lung cancer, with its low survival rates, is the biggest cancer killer in the UK – Over one fifth of all cancer deaths in UK are from lung cancer (25% for male and 18% for female) • Major risk factor is smoking. – Smoking time trends different for men/women: uptake of smoking started to decrease in cohorts of men after 1970, while for women the levelling off was later, after 1980 • Other risk factors include exposure to workplace agents, radon, air pollution … Interested in similarity and specificity of patterns between men and women Space-time analysis of Male and Female lung cancer incidence Male/Female lung cancer incidence in Yorkshire:81-85, 86-90, 91-95, 96-99 (Richardson, Abellan, Best, 2006) Shared and specific patterns and time trends Shared component Female/Male differential Time trend for male RRs in 10 wards Time trend for female RRs in 10 wards Detection of space-time interaction patterns Detection of space-time interaction patterns Noisy data in each area Noise model: Poisson/Binomial Latent structure: Space + Time + (Residuals) + joint Bayesian estimation Inference Detection of space-time interaction patterns Noisy data in each area Noise model: Poisson/Binomial Latent structure: Space + Time + Interactions Any patterns? + joint Bayesian estimation + Inference Detection of space-time interaction patterns • Study the persistence of patterns over time – Interpreted as associated with stable risk factors, environmental effects, socio-economic determinants • Highlight unusual patterns, via the inclusion of space time interaction terms, which are modelled by a mixture model • Unusual patterns in some areas may be linked to recording changes, emerging environmental hazards, impact of new policy or intervention program, … a general tool for surveillance ? Detection of space-time interaction patterns Yit ~ Poisson(Rit Eit); log Rit = Si + Tt + Uit Si ~ spatial CAR model (common spatial pattern) Tt ~ random walk (RW) model (common temporal trend) Uit ~ Normal(0, v) (space-time interaction; idiosyncratic variation) © Imperial College London Detection of space-time interaction patterns Yit ~ Poisson(Rit Eit); log Rit = Si + Tt + Uit Si ~ spatial CAR model (common spatial pattern) Tt ~ random walk (RW) model (common temporal trend) Uit ~ q Normal(0, v1) +(1-q) Normal(0, v2); v2 > v1 (mixture model to characterise ‘stable’ and ‘unstable’ patterns over time) • Compute posterior probability, pit, that interaction parameter Uit comes from the Normal(0, v2) component • Classify area as ‘unstable’ if pit > 0.5 for at least one time, t (simulation study → 10% false positive rate; 20% false negative rate) © Imperial College London Detecting unusual trends in congenital anomalies rates in England (Abellan et al 2008) • Annual postcoded data on congenital anomalies (non chromosomal) recorded in England for the period 1983 – 1998 • Annual postcoded data on total number of live births, still births and terminations • 136,000 congenital anomalies 84.5 per 105 birth-years • Congenital anomalies are sparse: Grid of 970 grid squares with variable size, to equalize the number of births and expected cases per area • Variations could be linked to socio-economic or environmental risk factors or heterogeneity in recording practices Interest in characterising space time patterns © Imperial College London Congenital anomalies in England, 1983-1998 Spatial main effect: evidence of spatial heterogeneity, linked to deprivation and maternal age Temporal main effect: downward trend around 1990 reflects implementation of “minor anomalies” exclusion policy Congenital anomalies: Space-time interactions Most areas are stable (cluster 1) Some have a change around 90-91 where modifications in the classification of anomalies occurred (clusters 2 and 3) Identified one very unusual time profile due to a change of local recording practice Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends Summary 6. Hierarchical related regression models for combining individual and small area data Summary • Smoothing of small area risks is important to help separate ‘signal’ (spatial pattern) from ‘noise’ –Allows meaningful inference even when data are sparse • Achieved by ‘borrowing’ information from neighbouring regions • Bayesian hierarchical modelling provides formal method for carrying out this ‘borrowing of information’ –Provides rich output for statistical inference (estimation, quantification of uncertainty, hypothesis testing) –But, depends on “structural” assumptions built into the model (e.g. spatial dependence) –Computationally intensive © Imperial College London Summary Bayesian approach extends naturally to allow: • Adjustment for covariates (see later) • Joint mapping of 2 or more health outcomes • Joint modelling of spatial and temporal variation Benefits of Joint Analysis of related health outcomes Joint analysis of two related health outcomes is of interest in several contexts: • Epidemiology: quantify ‘expected’ variability linked to shared risk factors and tease out specific patterns • Health planning: assess the performance of the health system, e.g. for health outcomes linked to screening policies • Data quality issues: uncover anomalous patterns linked to a data source shared by several outcomes Benefits of Space Time Analysis for (non-infectious) health outcomes • Study the persistence of patterns over time – Interpreted as associated with stable risk factors, environmental effects, distribution of health care access … • Highlight unusual patterns in time profiles via the inclusion of space-time interaction terms – Time localised excesses linked to e.g. emerging environmental hazards with short latency – Variability in recording practices Increased epidemiological interpretability Potential tool for surveillance Outline 1. Introduction 2. Mapping and spatial smoothing of health data 3. Classifying areas and summarising geographical variations in health outcomes 4. Modelling and mapping multiple outcomes 5. Modelling and mapping temporal trends 6. Hierarchical related regression models for combining individual and small area data Introduction • Models and applications discussed so far have focused on: – describing geographical and temporal patterns in health outcomes – partitioning sources of variation into e.g. systematic and idiosyncratic, spatial and temporal, shared and specific, … • Growing interest in trying to explain geographical variations at level of areas and individuals → Build regression models linking health outcomes and explanatory variables Regression models for small area data Standard regression for individual-level outcomes Individual exposure xij Aggregate exposure Zi, Xi yij Individual outcome Ecological regression Aggregate exposure Zi, Xi Yi Aggregate outcome Hierarchical Related Regression (HRR; Jackson et al, 2006, 2008a,b) Individual exposure Aggregate exposure xij Zi, Xi yij Yi Individual outcome Aggregate outcome Case study: Socioeconomic inequalities in health Jackson, Best and Richardson (2008b) • Geographical inequalities in health are well documented • One explanation is that people with similar characteristics cluster together, so area effects are just the result of differences in characteristics of people living in them (compositional effect) • But, evidence suggests that attributes of places may influence health over and above effects of individual risk factors (contextual effect) – economic, environmental, infrastructure, social capital/cohension Question • Is there evidence of contextual effects of area of residence on risk of limiting long term illness (LLTI) and heart disease, after adjusting for individual-level socio-demographic characteristics Data and Methodological Issues Methodological issues • Surveys typically contain sparse individual data per area so difficult to estimate contextual effects • Can’t separate individual and contextual effects using only aggregate data (ecological bias) • Improve power and reduce bias by combining data using new class of multilevel models developed by BIAS Our goal: data synthesis using • Individual-level data – Health Survey for England, 1997-2001 • Area-level (electoral ward) data – 1991 census small-area statistics – Hospital Episode Statistics Data sources INDIVIDUAL DATA Health Survey for England AREA (WARD) DATA Census small area statistics • Self-reported limiting long term illness • Self reported hospitalisation for heart disease • Carstairs deprivation index (area-level material deprivation) • age and sex • ethnicity • social class • car access • income • etc. Ward codes made available under special license Individual-level Health outcomes Contextual effect Individual predictors Multilevel model for individual data yij = disease (1) / no disease (0) xij = non-white (1) / white (0) Zi = deprivation score b c m,v2 ai xij yij person j Zi area i yij ~ Bernoulli(pij), person j, area i logit pij = ai + b xij + c Zi ai ~ Normal(m, v2) b = relative risk of disease for non-white versus white individual c = contextual effects ai = “unexplained” area effects Results from analysis of individual survey data: Heart Disease (n=5226) Univariate regression Area deprivation Multiple regression Area deprivation No car Social class IV/V Non white 0.2 0.5 1.0 2.0 odds ratio 5.0 10 Results from analysis of individual survey data: Limiting Long Term Illness (n=1155) Univariate regression Area deprivation Multiple regression Area deprivation Female Non white Doubled income 0.5 1.0 odds ratio 2 3 Comments • CI wide and not significant for most effects • Some evidence of contextual effect of area deprivation for both heart disease and LLTI – Adjusting for individual risk factors (compositional effects) appears to explain contextual effect for heart disease – Unclear whether contextual effect remains for LLTI after adjustment for individual factors – Survey data lack power to provide reliable answers about contextual effects • What can we learn from aggregate data? Area-level data AREA (WARD) DATA Contextual effect Census small area statistics • Carstairs deprivation index • population count by age and sex • proportion reporting LLTI Aggregate health outcomes Aggregate versions &of denominators individual predictors • proportion non-white • proportion in social class IV/V • proportion with no car access PayCheck (CACI) • mean & variance of household income Hospital Episode Statistics • number of admissions for heart disease Standard ecological regression model Yi = number with disease Ni = population Xi = proportion non-white Zi = area deprivation score Yi ~ Binomial(qi, Ni), area i logit qi = Ai + BXi + CZi Ai ~ Normal(M, V2) B = association between disease prevalence and proportion non-white C = contextual effects Ai = “unexplained” area effects B M,V2 C Ai Zi area i Yi Xi Ni Comparison of individual and ecological regressions: Heart Disease Individual Area deprivation Ecological No car Social class IV/V Non white 0.2 0.5 2.0 1.0 odds ratio 5.0 10 Comparison of individual and ecological regressions: Limiting Long Term Illness Individual Area deprivation Ecological Female Non white Doubled income 0.5 1.0 2 odds ratio 3 6 Ecological bias Ecological bias (difference between individual and aggregate level effects) can be caused by: • Confounding – confounders can be area-level (between-area) or individuallevel (within-area). → include control variables and/or random effects in model • Non-linear covariate-outcome relationship, combined with within-area variability of covariate – No bias if covariate is constant in area (contextual effect) – Bias increases as within-area variability increases – …unless models are refined to account for this hidden variability Standard ecological regression model Yi = number with disease Ni = population Xi = proportion non-white Zi = area deprivation score Yi ~ Binomial(qi, Ni), area i logit qi = Ai + BXi + CZi Ai ~ Normal(M, V2) B = association between disease prevalence and proportion non-white C = contextual effects Ai = “unexplained” area effects B M,V2 C Ai Zi area i Yi Xi Ni Integrated ecological regression model Yi = number with disease Ni = population Average of the Xi = proportionindividual non-whiteprobabilities of disease, Zi = area deprivation scorepij, in area i 2 m,v Yi ~ Binomial(qi, Ni), b c area i qi = pij(xij,Zi,ai, b, c)fi(x)dx ai ~ Normal(m, v2) b = relative risk of disease for non-white versus white individual c = contextual effects ai = “unexplained” area effects ai Zi area i Yi Xi Ni Combining individual and aggregate data Multilevel model for individual data b xij c yij person j Integrated ecological model m,v2 m,v2 ai ai Zi Zi area i area i b Yi c Xi Ni Combining individual and aggregate data Hierarchical Related Regression (HRR) model b c m,v2 Joint likelihood for yij and Yi depending on shared parameters ai, b, c, m, v2 (Jackson, Best, Richardson, 2006, 2008a,b) ai xij yij Zi person j area i Yi Xi Ni Combining individual and aggregate data Hierarchical Related Regression (HRR) model b c m,v2 Joint likelihood for yij and Yi depending on (Jackson, Best, Richardson, Estimation carried out using shared parameters 2006, 2008a,b) ai, b, c, m, v2 R software (maximum likelihood) or WinBUGS (Bayesian) ai xij yij Zi person j area i Yi Xi Ni Comparison of results from different regression models: Heart Disease Individual Area deprivation Standard ecological Integrated ecological No car HRR Social class IV/V PR95 = 10.1; 95% CI(5.3, 18.1) Non white PR95 = 4.2; 95% CI(3.6, 5.1) 0.2 0.5 1.0 2.0 odds ratio 5.0 10 Comparison of results from different regression models: Limiting Long Term Illness Individual Area deprivation Standard ecological Integrated ecological Female HRR PR75 = 2.7; 95% CI(1.7, 4.1) Non white PR75 = 2.9; 95% CI(2.4, 3.7) Doubled income 0.5 1.0 odds ratio 2 3 6 Comments • Integrated ecological model yields odds ratios that are consistent with individual level estimates from survey • Large gains in precision achieved by using aggregate data • Significant contextual effect of area deprivation for LLTI but not heart disease • More unexplained between-area variation (PR95) for heart disease than LLTI • Little difference between estimates based on aggregate data alone and combined individual + aggregate data – Individual sample size very small (~0.1% of population represented by aggregate data) • In other applications with larger individual sample sizes and/or less informative aggregate data, combined HRR model yields greater improvements (simulation study) Strengths of HRR approach…… Aims to provide individual-level inference using aggregate data by: Fitting integrated individual-level model to alleviate one source of ecological bias Including samples of individual data to help identify effects Uses data from all geographic areas (wards, constituencies), not just those in the survey Improves precision of parameter estimates Improves ability to investigate contextual effects …..and limitations of HRR approach Integrated individual-level model relies on large contrasts in the predictor proportions across areas e.g. limited variation in % non-white across constituencies: (median 2.7%, 95th percentile 33) Our estimates may not be completely free from ecological bias (Jackson et al, 2006) If individual level data too sparse, may be overwhelmed by aggregate data Data requirements for HRR models • Individual data – requires geographical (group) identifiers for individual data • Aggregate data – requires large exposure contrasts between areas – requires information on within-area distribution of covariates • Important to check compatibility of different data sources when combining data Thank you for your attention Acknowledgements: Sylvia Richardson, Juanjo Abellan, Virgilio Gomez-Rubio, Chris Jackson Training courses in Bayesian Analysis of Small Area Data using WinBUGS and INLA, London, July 13-16 2010 See www.bias-project.org.uk for details References Abellan JJ, Richardson S and Best N. Use of space-time models to investigate the stability of patterns of disease. Environ Health Perspect 116(8), (2008), 1111-1119. Best N and Hansel A. Geographic variations in risk: adjusting for unmeasured confounders through joint modelling of multiple diseases. Epidemiology, 20(3), (2009), 400-410. Best, N.G., Richardson, S. and Thomson, A. A comparison of Bayesian spatial models for disease mapping. Statistical Methods in Medical Research 14, (2005), 35-59. Jackson C, Best N and Richardson S. Hierarchical related regression for combining aggregate and individual data in studies of socio-economic disease risk factors. J Royal Statistical Society Series A: Statistics in Society 171 (2008a) ,159-178 Jackson C, Best N and Richardson S. Studying place effects on health by synthesising individual and arealevel outcomes. Social Science and Medicine, 67, (2008b), 1995-2006 Jackson C, Best N and Richardson S. Improving ecological inference using individual-level data. Statistics in Medicine, 25, (2006), 2136-2159 Knorr-Held, L. and Best, N.G., A shared component model for detecting joint and selective clustering of two diseases, Journal of the Royal Statistical Society, Series A 164, (2001), 73-86. Richardson, S., Thomson, A., Best, N.G. and Elliott, P. Interpreting posterior relative risk estimates in disease mapping studies. Environmental Health Perspectives 112, (2004), 1016-1025. Richardson, S., Abelan, J.J.,and Best, N. Bayseian spatio-temporal analysis of joint patterns of male and female lung cancer in Yorkshire (UK). Statistical Methods in Medical Research 15, (2006), 385-407. Tzala, E. Best, N. Bayesian latent variable modelling of multivariate spatio-temporal variation in cancer mortality. Stat Methods Med Res 17, (2008), 97-118. Spanish disease atlases and other related resources López-Abente G, Pollan M, Escolar A, Errezola M, Abraira V. 1996. Atlas of cancer mortality and other causes of death in Spain 1978-1992. Madrid: Fundación Científica de la Asociación Española contra el Cáncer. http://www2.uca.es/hospital/atlas92/www/Atlas92.html. Maps of SMRs and Directly Standardised Rates at province level Benach J, Yasui Y, Borrell C, Rosa E, Pasarin MI, Benach N, Español M, Martinez JM, Daponte A. 2001. Atlas of mortality at small area level in Spain 1987-1995. Barcelona: Universitat Pompeu Fabra. Empirical Bayes smoothing of relative risks at municipality level (aggregated some municipalities to reduced sparseness). All of the following produce maps of Bayesian smoothed relative risks at municipality level: López-Abente G, Ramis R, Pollán M, Aragonés N, Pérez-Gómez B, Gómez-Barroso D, Carrasco JM, Lope V, García-Pérez J, Boldo E, García-Mendizabal MJ. 2007. Atlas of cancer mortality at municipality level in Spain 1989-1998. Área de Epidemiología Ambiental y Cáncer del Centro Nacional de Epidemiología, ISCIII. Botella P, Zurriaga O, Posada de la Paz M, Martinez-Beneito MA, Bel E, Robustillo A, Ramalle E, Duran E, Sanchez-Porro P. 2006. National-Provincial Atlas of Rare Diseases 1999-2003. Martinez-Beneito MA, Lopez-Quilez A, Amador A, Melchor I, Botella P, Abellan C, Abellan JJ, Verdejo F, Zurriaga O, Vanaclocha H, Escolano M. 2005. Atlas of Mortality in the Valencian Region 1991-2000. Benach J, Matinez JM, Yasui Y, Borrell C, Pasarin MI, Español E, Benach N. 2004. Atlas of mortality at small area level in Catalonia 1984-1998. Barcelona: Universitat Pompeu Fabra / Fundació Jaume Bofill / Editorial Mediterrània Spanish disease atlases and other related resources DEMAP group, Andalusian School of Public Health, Granada. Produced interactive mortality atlas for Andalucia, and socioeconomic indices at municipality level for Spain. See http://www.demap.es/Demap/index.html MEDEA: Research network on Epidemiology and Public Health, working on socioeconomic and environmental inequalities in health at small area level. See http://www.proyectomedea.org/medea.html VPM Atlas Project (Atlas de Variaciones en la Práctica Médica). Studying and mapping variations in provision and usage of health care at small area level in 16 of the 17 regions of Spain, using data on hospital discharges. See http://www.atlasvpm.org/avpm/inicio.inicio.do Smoothing of the RRs of hot spots (4 contiguous areas with average expected counts ≈ 5) for different spatial models Richardson et al (EHP, 2004) True RR = 3 True RR = 2