1. No category

Analysis of time-stratified casecrossover studies in environmental
epidemiology using Stata
Aurelio Tobías
Spanish Council for Scientific Research (CSIC), Barcelona, Spain
Ben Armstrong and Antonio Gasparrini
London School of Hygiene and Tropical Medicine, UK
2014 UK Stata Users Group meeting
London, 12th September 2014
Background
• The time stratified case crossover design is a
popular alternative to conventional time
series regression for analysing associations
between time series of environmental
exposures (air pollution, weather) and counts
of health outcomes
The case-crossover design
• Proposed by Maclure (1991) to study transient
effects on the risk of acute health events
… compared to your usual
routine?
Have you
made ​any
unusual activity
…
Health
event
The case-crossover design
• Proposed by Maclure (1991) to study transient
effects on the risk of acute health events
Exposed?
Control period
Exposed?
Risk period
Health
event
The case-crossover design
• Proposed by Maclure (1991) to study transient
effects on the risk of acute health events
Exposed?
Control period
• Analysis likewise a
matched case-control
Exposed?
Risk period
Health
event
Control
Exp. No
Risk Exp.
period No
a
c
b
d
OR = b/c
Application to environmental studies
• Firstly adapted to air pollution studies
– Philadelphia (Neas et al. 1999) and Barcelona (Sunyer et al.
2000)
• comparing exposure levels for a given day (t) when
health event occurs vs. levels before (t-7) and after
(t+7) the health event
• Allows to control for time-trend and seasonality by
design, since compares exposure levels between
same weekdays within each month of each year
Time-stratified case-crossover
(Figueiras et al. 2010)
Option 1: conditional logistic regression
• Dataset needs to be reshaped from time-series
format to individual matched case-control by using
the new command,
• . mkcco vary, date(vardate)
–
–
–
–
This also creates three new variables,
case case days (=1)
wn
weight observations by n events on case day
ccset set num. to match case-control days
• Then, use Stata’s clogit command
. use madrid, clear
. list date y pm10 temp in 1/31, noobs clean
date
01jan2003
02jan2003
03jan2003
04jan2003
05jan2003
06jan2003
07jan2003
08jan2003
09jan2003
10jan2003
11jan2003
12jan2003
13jan2003
14jan2003
15jan2003
16jan2003
17jan2003
18jan2003
19jan2003
20jan2003
21jan2003
22jan2003
23jan2003
24jan2003
25jan2003
26jan2003
27jan2003
28jan2003
29jan2003
30jan2003
31jan2003
y
64
56
73
69
71
68
63
85
67
80
65
95
60
76
77
75
76
75
74
79
73
72
72
74
59
77
64
65
74
77
55
pm10
14
12
16
14
17
9
19
13
13
22
23
17
45
60
72
54
65
43
12
19
16
21
25
46
42
26
22
41
18
17
17
temp
9.6
11
9.8
8.8
4.3
7
3
6.8
4.6
2.7
1
.85
.5
.45
1.3
2.6
2
.8
6.3
5.8
8
6.7
7
6.2
8.1
10
15
9.7
5.7
4.8
2.3
Sun
Mon
Tue
Wed
Thu
Fri
Sat
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
. generate dow = dow(date)
. list date y pm10 temp in 1/31 if dow==1, noobs clean
date
06jan2003
13jan2003
20jan2003
27jan2003
y
68
60
79
64
pm10
9
45
19
22
temp
7
.5
5.8
15
. mkcco y, date(date)
Dataset reshaped from 1096 to 5480 obs. (1096 case days and 4384 control days).
New variables case, wn and ccset have been added to the resaphed dataset.
Sun
Mon
Tue
Wed
Thu
Fri
Sat
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
. list date y pm10 temp case ccset wn in 1/36 if dow==1, noobs clean
date
06jan2003
13jan2003
20jan2003
27jan2003
...
y
68
60
79
64
pm10
9
45
19
22
temp
7
.5
5.8
14.8
case
1
0
0
0
ccset
2003021
2003021
2003021
2003021
wn
68
68
68
68
Sun
Mon
Tue
Wed
Thu
Fri
Sat
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
...
date
06jan2003
13jan2003
20jan2003
27jan2003
y
68
60
79
64
pm10
9
45
19
22
temp
7
.5
5.8
14.8
case
0
1
0
0
ccset
2003022
2003022
2003022
2003022
wn
60
60
60
60
Sun
Mon
Tue
Wed
Thu
Fri
Sat
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
...
date
06jan2003
13jan2003
20jan2003
27jan2003
y
68
60
79
64
pm10
9
45
19
22
temp
7
.5
5.8
14.8
case
0
0
1
0
ccset
2003023
2003023
2003023
2003023
wn
79
79
79
79
Sun
Mon
Tue
Wed
Thu
Fri
Sat
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
...
date
06jan2003
13jan2003
20jan2003
27jan2003
y
68
60
79
64
pm10
9
45
19
22
temp
7
.5
5.8
14.8
case
0
0
0
1
ccset
2003024
2003024
2003024
2003024
wn
64
64
64
64
. clogit case pm10 [w=wn], group(ccset)
(frequency weights assumed)
...
Conditional (fixed-effects) logistic regression
Log likelihood =
-98913.41
Number of obs
LR chi2(1)
Prob > chi2
Pseudo R2
=
=
=
=
295025
8.39
0.0038
0.0000
-----------------------------------------------------------------------------case |
Coef.
Std. Err.
z
P>|z|
[95% Conf. Interval]
-------------+---------------------------------------------------------------pm10 |
.0007406
.0002552
2.90
0.004
.0002404
.0012408
------------------------------------------------------------------------------
Option 1: conditional logistic regression
• Advantage
– Similar approach to a matched
case-control study, being
familiar to the epidemiologist
• Limitations
– Unfriendly data management
– Large (huge) data-sets at
individual level
– Slow convergence when fitting
interaction terms for
individual factors (sex, age)
– Cannot account for
overdispersion and residual
autocorrelation
Option 2: Poisson regression
• Lu and Zeger (2007)
– “ … the time-stratified case-crossover design leads to the
same estimate as obtained from a Poisson regression with
dummy variables indicating the strata”
• Strata defined by the 3-way interaction between
year, month and day of week
• Either poisson or glm (jointly with the scale
option) Stata’s commands can be used
.
.
.
.
use madrid, clear
generate yy = year(date)
generate mm = month(date)
generate dow = dow(date)
. poisson y pm10 i.yy##i.mm#i.dow
...
Poisson regression
Log likelihood =
-3787.589
Number of obs
LR chi2(252)
Prob > chi2
Pseudo R2
=
=
=
=
1096
1374.12
0.0000
0.1535
-----------------------------------------------------------------------------y |
Coef.
Std. Err.
z
P>|z|
[95% Conf. Interval]
-------------+---------------------------------------------------------------pm10 |
.0007406
.0002552
2.90
0.004
.0002404
.0012408
|
...
252 parameters!
|
_cons |
4.35927
.0563535
77.36
0.000
4.248819
4.469721
------------------------------------------------------------------------------
. glm y pm10 i.yy##i.mm#i.dow, family(poisson) link(log) scale(x2)
...
Generalized linear models
Optimization
: ML
Deviance
Pearson
=
=
1069.73306
1065.616046
Variance function: V(u) = u
Link function
: g(u) = ln(u)
No. of obs
Residual df
Scale parameter
(1/df) Deviance
(1/df) Pearson
=
=
=
=
=
1096
843
1
1.26896
1.264076
[Poisson]
[Log]
AIC
= 7.373338
Log likelihood
= -3787.588997
BIC
= -4830.78
-----------------------------------------------------------------------------|
OIM
y |
Coef.
Std. Err.
z
P>|z|
[95% Conf. Interval]
-------------+---------------------------------------------------------------pm10 |
.0007406
.0002869
2.58
0.010
.0001782
.001303
|
...
|
_cons |
4.35927
.0633589
68.80
0.000
4.235088
4.483451
-----------------------------------------------------------------------------(Standard errors scaled using square root of Pearson X2-based dispersion.)
Option 2: Poisson regression
• Advantages
– Avoids reshape of dataset
keeping original time-series
format
– Can account for
overdispersion and residual
autocorrelation
• Limitations
– Large number of interaction
parameters, slowing the
convergence time (even
making difficult it with missing
data)
– Even slower when fitting
interaction terms for
individual factors (sex, age)
Option 3: Conditional Poisson regression
• Armstrong and Gasparrini (ISEE 2011)
– “ … Poisson models with large numbers of indicator variables
can alternatively be fit with conditional Poisson models,
conditioning on numbers of events in the time stratum”
– Firstly proposed by Farrington (1995) for the self-controlled
case-series design for vaccine safety studies
• Strata defined by grouping same weekday within each
month of each year, and then use Stata’s xtpoisson
command jointly with the fe option
• But overdispersion needs to accounted for with the new
xtpois_ovdp command
. egen set = group(yy mm dow)
. xtset set date
...
. xtpoisson y pm10, fe
...
Conditional fixed-effects Poisson regression
Group variable: set
Log likelihood
= -2854.4979
Number of obs
Number of groups
=
=
1096
252
Obs per group: min =
avg =
max =
4
4.3
5
Wald chi2(1)
Prob > chi2
=
=
8.42
0.0037
-----------------------------------------------------------------------------y |
Coef.
Std. Err.
z
P>|z|
[95% Conf. Interval]
-------------+---------------------------------------------------------------pm10 |
.0007406
.0002552
2.90
0.004
.0002404
.0012408
-----------------------------------------------------------------------------. xtpois_ovdp
Estimate and standard errors corrected for over-dipersion
df: 843 ; pearson x2: 1065.6 ; dispersion:
1.26
-----------------------------------------------------------------------------y |
Coef.
Std. Err.
z
P>|z|
[95% Conf. Interval]
-------------+---------------------------------------------------------------pm10 |
.0007406
.0002869
2.58
0.010
.0001782
.001303
------------------------------------------------------------------------------
CPU time
• My example (3 years data, 60 events/day)
– clogit
– poisson
– xtpoisson
0.181 sec.
2.111 sec.
0.066 sec.
Option 3: Conditional Poisson regression
• Advantages
– Avoids reshape of dataset
keeping original time-series
format
– Can account for
overdispersion and residual
autocorrelation
– Easy fit of interaction terms
for individual factors (sex, age)
with faster convergence
• Limitations
–?
Summary
• Case-crossover studies can easily be analysed using
Stata
1) Conditional logistic regression requires unfriendly data
management, and large (huge) data-sets at individual
level
2) Poisson, 3-way interaction, regression requires to fit large
number of interaction parameters, making it difficult
convergence
3) Conditional Poisson seems to solve both
Thanks
for your
attention!
Thanks
for your
attention!