Analysing a 2x2 Crossover: Separating Treatment From Period and Carry-Over

On this page
อ่านฉบับภาษาไทย (Thai version)
Abstract
A 2x2 crossover trial gives every participant both treatments, one per period, in order AB or BA. Readings move with treatment, with a period effect shared by everyone in period 2, and possibly with carry-over: the first drug still acting after the switch. Half the difference between the two sequences' mean period differences (each participant's period-1 reading minus period-2 reading) estimates the treatment effect free of any period effect. Carry-over is inseparable from the sequence comparison. In a simulated trial of 24 adults with hypertension, the estimated effect of drug A versus drug B on systolic blood pressure is -8.02 mmHg (95% CI -11.74 to -4.30). With the simulated carry-over of -2 mmHg the estimator targets -7 mmHg; this draw lies near the true -8 by chance. The carry-over test has a power (chance of detecting a real effect) of 0.054. A washout, a treatment-free gap long enough for the first drug to fade, is the defence against carry-over, not a preliminary test.
Two orders, two periods, one drug comparison
Twenty-four adults with hypertension join a trial of two fictional blood pressure drugs, A and B. Each takes both drugs in turn, with systolic blood pressure measured at the end of each treatment period. Twelve are randomised to the order AB and twelve to the order BA.
At the steering meeting, the interim summary shows the two orders moving in opposite directions after the switch (simulated data). Mean systolic blood pressure rises by 3.88 mmHg from period 1 to period 2 in the group that took A first (sequence AB) and falls by 12.15 mmHg in the group that took B first (sequence BA). Averaged over both orders, period 2 sits lower.
Part of each change belongs to the drug, part may belong to period 2 itself, and part to whatever the first drug left behind. The committee wants one number for the effect of A versus B. The analysis must separate that drug effect from the period and from any effect of the first drug that lingers after the switch.
What a crossover trial is, and what moves its readings
A crossover trial gives each participant every treatment under study, one at a time, in a randomised order. Each treatment interval is a period, and the order a participant follows is that participant's sequence. The 2x2 design has two periods and two sequences, AB and BA.
The design rests on a within-participant comparison: each person is measured on both drugs, so stable differences between people, such as their usual blood pressure, cancel. When those differences are large, a crossover can match the precision of a much larger parallel-group trial, in which each participant receives only one treatment [1, 2]. The price is that time and treatment order now enter the data.
Four effects can be named. The treatment effect $\Delta$ (delta) is the difference between drug A and drug B within a participant. The period effect $\pi$ (pi) is a shift shared by everyone in period 2, whichever drug they take.
The carry-over effect $\lambda$ (lambda) is an effect of the period-1 drug that persists into period 2. A washout period, a treatment-free interval between periods long enough for the first drug's effect to fade, is the design defence against carry-over.
A sequence effect is any systematic difference between the AB and BA groups as wholes. The table shows what each effect follows, and the next sections show why the last two cannot be separated.
| Effect | It follows | Symbol | Can it be estimated? |
|---|---|---|---|
| Treatment | the drug taken in the current period | $\Delta$ | Yes, within participants (shifted by minus half the carry-over, if any) |
| Period | whether the reading is in period 1 or period 2 | $\pi$ | Yes, within participants (shifted by half the carry-over, if any) |
| Carry-over | the drug taken in the previous period | $\lambda$ | Only between participants, and only imprecisely |
| Sequence | the order, AB or BA | none of its own | No: the design gives it the same pattern in the data as carry-over (the two are aliased), so no analysis can separate them |
A model for the four cell means
Write $Y_{ij}$ for the outcome of participant $i$ in period $j$, where $j$ is 1 or 2. The design has four cells, one for each sequence in each period, and a cell mean is the average reading in one cell. The standard model for a 2x2 crossover is
$$Y_{ij} = \mu + s_i + \pi\,[\text{period 2}] + \Delta\,[\text{on A}] + \lambda\,[\text{period 2 of AB}] + e_{ij}$$Each bracket equals 1 when its condition holds and 0 otherwise. The constant $\mu$ (mu) is the mean on drug B in period 1, $s_i$ is participant $i$'s own level relative to it, and $e_{ij}$ is the noise of a single reading. Treating $s_i$ as a draw from a normal distribution makes it a random intercept, the term a mixed model (a regression with a random level for each participant) uses for each participant's own level.
Strictly, $\lambda$ is the carry-over of A minus the carry-over of B. Carry-over common to both drugs looks exactly like a period effect, so the model carries only the difference, placed on period 2 of sequence AB.
| Sequence | Period 1 | Period 2 |
|---|---|---|
| AB | $\mu + \Delta$ (on A) | $\mu + \pi + \lambda$ (on B) |
| BA | $\mu$ (on B) | $\mu + \pi + \Delta$ (on A) |
Period differences: the estimator
For each participant, subtract the period-2 reading from the period-1 reading: $d_i = Y_{i1} - Y_{i2}$. The participant's own level $s_i$ appears in both readings and cancels. Write $\bar d_{AB}$ and $\bar d_{BA}$ for the mean of $d_i$ in each sequence.
From the cell means, $\bar d_{AB}$ estimates $\Delta - \pi - \lambda$ and $\bar d_{BA}$ estimates $-\Delta - \pi$. Their difference removes $\pi$, and their sum removes $\Delta$:
$$\hat\Delta = \tfrac{1}{2}\left(\bar d_{AB} - \bar d_{BA}\right)$$The treatment estimate $\hat\Delta$ (delta-hat; a hat marks an estimate) is half the difference between the two sequences' mean period differences.
$$\hat\pi = -\tfrac{1}{2}\left(\bar d_{AB} + \bar d_{BA}\right)$$The period estimate $\hat\pi$ is minus half their sum, the shift from period 1 to period 2. With no carry-over ($\lambda = 0$), $\hat\Delta$ estimates $\Delta$ and $\hat\pi$ estimates $\pi$.
The standard error (SE) of $\hat\Delta$ is half the standard error of a two-sample t-test comparing $d_i$ between the sequences (Hills and Armitage [3]). The SE is the spread of the estimate over many repetitions of the same trial. That t-test estimates $\bar d_{AB} - \bar d_{BA}$, which is twice $\hat\Delta$, so its t statistic and P value apply to $\hat\Delta$ unchanged.
Hand example: a period effect with no drug difference
This hand example uses a symptom score from 0 to 10, where higher is worse. In sequence AB the mean score is 8 on A in period 1 and 5 on B in period 2. In sequence BA it is 8 on B in period 1 and 5 on A in period 2.
-
Mean period difference in AB
\[ \bar d_{AB} = 8 - 5 = 3 \]
Period 1 minus period 2 in sequence AB.
-
Mean period difference in BA
\[ \bar d_{BA} = 8 - 5 = 3 \]
The same subtraction in sequence BA.
-
Treatment effect
\[ \hat\Delta = \tfrac{1}{2}(3 - 3) = 0 \]
Drug A and drug B do not differ.
-
Period effect
\[ \hat\pi = -\tfrac{1}{2}(3 + 3) = -3 \]
Scores are 3 points lower in period 2, whichever drug is taken.
Result: The whole fall from 8 to 5 belongs to period 2, and none of it to either drug.
Inside sequence AB alone, the fall coincides with the switch to B and would be credited to B. Inside sequence BA alone, the same fall would be credited to A. Only the two sequences together separate the drug from the period.
Carry-over is aliased with sequence
Now let the first drug linger. Put the cell means into the two estimators and take expectations, written $E(\cdot)$ for the average over many repetitions of the same trial:
$$E(\hat\Delta) = \Delta - \tfrac{\lambda}{2}$$So the treatment estimate is shifted by minus half the carry-over. The period estimate is shifted the other way, $E(\hat\pi) = \pi + \lambda/2$.
The two mean period differences are the only within-participant information, and they give two equations, $\Delta - \pi - \lambda$ and $-\Delta - \pi$, in three unknowns. No combination of them isolates $\Delta$, or $\pi$, free of $\lambda$. The one estimate of $\Delta$ that is free of $\lambda$ compares period 1 alone (AB period 1 minus BA period 1), a between-participant comparison; the model with a carry-over term, later in this article, does exactly that.
Two effects are aliased when the design gives them the same pattern in the data, so no analysis can separate them.
In the 2x2 crossover, carry-over, a sequence effect and a treatment-by-period interaction (a drug difference that changes between periods) all show up as one quantity: the difference between the sequences in their participant totals (each participant's period-1 reading plus period-2 reading) [1, 2]. The data cannot say which of the three it measures.
The simulated trial below was generated with $\Delta = -8$, $\pi = -2$ and $\lambda = -2$ mmHg. The treatment estimate therefore targets $-8 - (-2)/2 = -7$ mmHg, a bias of 1 mmHg towards no effect. The period estimate targets $-2 + (-2)/2 = -3$ mmHg.
Why the carry-over test has little power
Carry-over can only be estimated from the participant totals, $t_i = Y_{i1} + Y_{i2}$. Their means differ between the sequences by $\lambda$, and a two-sample t-test on the totals is the carry-over test. Each total, however, contains $2s_i$, so the whole between-participant variation stays in it:
$$\operatorname{SE}(\hat\lambda) = \sqrt{\left(4\sigma_s^2 + 2\sigma_e^2\right)\left(\tfrac{1}{n_{AB}} + \tfrac{1}{n_{BA}}\right)}$$Here $\sigma_s$ is the between-participant standard deviation (SD), the spread of the levels $s_i$, and $\sigma_e$ is the within-participant SD of a single reading. The counts $n_{AB}$ and $n_{BA}$ are the participants in each sequence. The treatment estimate works on differences, where $\sigma_s$ cancels, so its standard error depends on $\sigma_e$ alone.
With the design values of the simulated trial, $\sigma_s = 12$ mmHg, $\sigma_e = 6$ mmHg and 12 participants per sequence, the standard error of the carry-over estimate is 10.39 mmHg.
Power is the probability that a test declares an effect when an effect of a stated size truly exists. The power of a two-sided test at the 5% level against the true $\lambda = -2$ mmHg is 0.054, barely above the 5% it would have with no carry-over at all. The carry-over that this design would detect in four trials out of five, a usual planning target for power, is about 30.5 mmHg, far larger than the treatment effect.
The carry-over test compares participant totals between sequences, a between-participant comparison with little power. A non-significant result does not show that carry-over is absent; carry-over is prevented by an adequate washout, not removed by a preliminary test.
Testing first, then using period 1 only: why the two-stage procedure fails
A once-standard approach, proposed by Grizzle [4], runs in two stages. First, carry-over is tested on the participant totals at a lenient significance level. If that test is significant, period 2 is discarded and the drugs are compared in period 1 as a parallel-group trial; otherwise, the period differences are used.
Freeman [5] showed that this procedure inflates the type I error, the probability of declaring a drug difference when none exists, and biases the final estimate. The carry-over test and the period-1 comparison share the period-1 readings and are strongly correlated, so the first stage selects which result the second stage reports.
Methodological work now advises against the two-stage procedure [1, 5]. Carry-over is handled at the design stage instead: a condition that is stable over the trial, a washout period between treatments long enough for the first drug's effect to fade, and outcomes measured at the end of each period. The crossover design and analysis guide covers the design side, including the washout, in more depth.
The analysis model: a mixed model with treatment and period
The usual analysis is a linear mixed model: a regression with fixed terms for treatment and period and a random intercept for each participant. It is fitted by REML (restricted maximum likelihood, the usual way to estimate the spread between and within participants in small samples).
The degrees of freedom of a t-test set the exact shape of the t distribution behind its P value and CI; Kenward-Roger degrees of freedom, a small-sample adjustment, are used for this model's t-tests.
With complete data, the treatment coefficient of this model equals $\hat\Delta$ exactly, and its standard error is half the standard error of the pooled-variance (Student) two-sample t-test on the period differences, the version that assumes one variance shared by both sequences. The t statistic, the P value and the 22 degrees of freedom are identical. A Welch t-test, which does not pool the two variances, would give slightly different degrees of freedom.
The model earns its place with less tidy data: it takes baseline covariates and keeps a participant with one missing period. The mixed-model series covers the model in depth.
Adding a carry-over term does not rescue the analysis. With $\lambda$ in the model, the treatment coefficient is estimated from period 1 alone, a between-participant comparison that gives up the crossover's main advantage.
The simulated trial in Stata and R
The dataset (simulated data) has 48 rows, one per participant per period: 24 participants, 12 in each sequence, with systolic blood pressure in mmHg. It was generated with $\Delta = -8$, $\pi = -2$ and $\lambda = -2$ mmHg, a between-participant SD of 12 mmHg and a within-participant SD of 6 mmHg.
The mean period differences are -3.88 mmHg in sequence AB and 12.15 mmHg in sequence BA. Half their difference gives $\hat\Delta = -8.02$ mmHg (standard error 1.79; 95% CI -11.74 to -4.30; P = 0.0002). The mixed model returns the same coefficient, standard error and interval. The period estimate is -4.13 mmHg (95% CI -7.85 to -0.41).
In this one draw, $\hat\Delta$ happens to lie near the true $\Delta$ of -8, although its target is -7. A single trial of 24 participants cannot reveal a bias of 1 mmHg against a standard error of 1.79.
The carry-over estimate is 7.45 mmHg (standard error 5.45; 95% CI -3.85 to 18.75; P = 0.19), against a true value of -2 mmHg. It is far from the truth and not significant, exactly what low power predicts. Its standard error is about half the design value of 10.39 because the fitted between-participant SD in this draw, 5.19 mmHg, is well below the 12 mmHg used to generate the data.
The model with a carry-over term estimates the treatment effect from period 1 only: -4.29 mmHg (95% CI -10.89 to 2.31).
| Quantity | Estimate (mmHg) | 95% CI (mmHg) | Target (mmHg) |
|---|---|---|---|
| Mean period difference, sequence AB | -3.88 | not shown | -4 |
| Mean period difference, sequence BA | 12.15 | not shown | 10 |
| Treatment effect, period differences | -8.02 | -11.74 to -4.30 | -7 |
| Treatment effect, mixed model | -8.02 | -11.74 to -4.30 | -7 |
| Period effect | -4.13 | -7.85 to -0.41 | -3 |
| Carry-over, totals AB minus BA | 7.45 | -3.85 to 18.75 | -2 |
| Treatment effect, model with carry-over term | -4.29 | -10.89 to 2.31 | -8 |
Stata: the whole analysis, from the cell means to the power of the carry-over test
* Analysis of a 2x2 AB/BA crossover (simulated data, not evidence about any real drug)
* Outcome sbp = systolic blood pressure in mmHg; drug A versus drug B; 12 + 12 participants, two periods.
version 18
clear all
set more off
set linesize 120
* read the simulated dataset (not published with this article)
import delimited using crossover.csv, clear varnames(1)
describe, short
list in 1/6, noobs sepby(id)
* participants per sequence (one row per participant)
egen byte first = tag(id)
tabulate seq if first
* cell means by sequence and period
table seq period, statistic(mean sbp) nformat(%9.2f)
* mixed model: treatment and period as fixed effects, a random intercept per participant
* (REML, Kenward-Roger degrees of freedom). Fitted quietly; the lines below print its treatment
* and period rows and the two standard deviations. R prints the same model in full.
quietly mixed sbp i.trt i.period || id:, reml dfmethod(kroger)
matrix T = r(table)
scalar mm_d = T[rownumb(T, "b"), colnumb(T, "sbp:1.trt")]
scalar mm_se = T[rownumb(T, "se"), colnumb(T, "sbp:1.trt")]
scalar mm_lo = T[rownumb(T, "ll"), colnumb(T, "sbp:1.trt")]
scalar mm_hi = T[rownumb(T, "ul"), colnumb(T, "sbp:1.trt")]
scalar mm_df = T[rownumb(T, "df"), colnumb(T, "sbp:1.trt")]
scalar mm_p = T[rownumb(T, "pvalue"), colnumb(T, "sbp:1.trt")]
scalar mm_pi = T[rownumb(T, "b"), colnumb(T, "sbp:2.period")]
scalar mm_pilo = T[rownumb(T, "ll"), colnumb(T, "sbp:2.period")]
scalar mm_pihi = T[rownumb(T, "ul"), colnumb(T, "sbp:2.period")]
scalar mm_sdid = exp(_b[lns1_1_1:_cons])
scalar mm_sdres = exp(_b[lnsig_e:_cons])
display "mixed model, treatment (A minus B): " %8.4f mm_d " SE " %7.4f mm_se " 95% CI " %8.4f mm_lo " to " %8.4f mm_hi
display "mixed model, treatment: df " %7.4f mm_df " P " %6.4f mm_p
display "mixed model, period (2 minus 1): " %8.4f mm_pi " 95% CI " %8.4f mm_pilo " to " %8.4f mm_pihi
display "mixed model, SD between participants " %7.4f mm_sdid " SD within participants " %7.4f mm_sdres
* mixed model with an explicit carry-over term: treatment is then estimated from period 1 only
quietly mixed sbp i.trt i.period i.carry || id:, reml dfmethod(kroger)
matrix C = r(table)
scalar mmc_d = C[rownumb(C, "b"), colnumb(C, "sbp:1.trt")]
scalar mmc_lo = C[rownumb(C, "ll"), colnumb(C, "sbp:1.trt")]
scalar mmc_hi = C[rownumb(C, "ul"), colnumb(C, "sbp:1.trt")]
scalar mmc_l = C[rownumb(C, "b"), colnumb(C, "sbp:1.carry")]
scalar mmc_lse = C[rownumb(C, "se"), colnumb(C, "sbp:1.carry")]
scalar mmc_lp = C[rownumb(C, "pvalue"), colnumb(C, "sbp:1.carry")]
display "carry-over model, treatment (A minus B): " %8.4f mmc_d " 95% CI " %8.4f mmc_lo " to " %8.4f mmc_hi
display "carry-over model, carry-over term: " %8.4f mmc_l " SE " %7.4f mmc_lse " P " %6.4f mmc_lp
* one row per participant: period difference d = period 1 minus period 2, and the participant total
keep id seq period sbp
reshape wide sbp, i(id) j(period)
generate double d = sbp1 - sbp2
generate double tot = sbp1 + sbp2
generate byte grp = cond(seq == "AB", 1, 2)
label define grp 1 "AB" 2 "BA"
label values grp grp
* treatment effect from period differences: Delta-hat = (mean d in AB - mean d in BA)/2, two-sample t-test
ttest d, by(grp)
scalar dbar_ab = r(mu_1)
scalar dbar_ba = r(mu_2)
scalar df_d = r(df_t)
scalar delta = (r(mu_1) - r(mu_2)) / 2
scalar delta_se = r(se) / 2
scalar delta_p = r(p)
scalar tc = invttail(df_d, 0.025)
display "treatment (A minus B): " %8.4f delta " SE " %7.4f delta_se " 95% CI " %8.4f delta - tc * delta_se " to " %8.4f delta + tc * delta_se " P " %6.4f delta_p
* period effect (period 2 minus period 1) from the same differences: -(mean d in AB + mean d in BA)/2
scalar per = -(dbar_ab + dbar_ba) / 2
display "period (2 minus 1): " %8.4f per " 95% CI " %8.4f per - tc * delta_se " to " %8.4f per + tc * delta_se
* carry-over test: compare participant totals between sequences (between-participant, so low power)
ttest tot, by(grp)
* design values used to simulate: sigma_s = 12, sigma_e = 6, lambda = -2
* power of the carry-over test at lambda = -2, and the carry-over it detects with 80 percent power
* (two-sided 5 percent level); SE of lambda-hat = sqrt((4 sigma_s^2 + 2 sigma_e^2) (1/n_AB + 1/n_BA))
scalar sig_s = 12
scalar sig_e = 6
scalar lam = -2
count if grp == 1
scalar n_ab = r(N)
count if grp == 2
scalar n_ba = r(N)
scalar des_se = sqrt((4 * sig_s^2 + 2 * sig_e^2) * (1 / n_ab + 1 / n_ba))
scalar des_df = n_ab + n_ba - 2
scalar des_t = invttail(des_df, 0.025)
scalar des_power = nt(des_df, lam / des_se, -des_t) + nttail(des_df, lam / des_se, des_t)
scalar des_mde80 = (des_t + invttail(des_df, 0.20)) * des_se
display "design SE of the carry-over estimate: " %7.4f des_se
display "power of the carry-over test at lambda = -2: " %6.4f des_power
display "carry-over detected with 80% power: " %7.4f des_mde80
. * Analysis of a 2x2 AB/BA crossover (simulated data, not evidence about any r
> eal drug)
. * Outcome sbp = systolic blood pressure in mmHg; drug A versus drug B; 12 + 1
> 2 participants, two periods.
. version 18
. clear all
. set more off
. set linesize 120
.
. * read the simulated dataset (not published with this article)
. import delimited using crossover.csv, clear varnames(1)
(encoding automatically selected: ISO-8859-1)
(7 vars, 48 obs)
. describe, short
Contains data
Observations: 48
Variables: 7
Sorted by:
Note: Dataset has changed since last saved.
. list in 1/6, noobs sepby(id)
+----------------------------------------------------+
| id seq period treatm~t trt carry sbp |
|----------------------------------------------------|
| 1 BA 1 B 0 0 146.2 |
| 1 BA 2 A 1 0 134.5 |
|----------------------------------------------------|
| 2 AB 1 A 1 0 152.4 |
| 2 AB 2 B 0 1 149.6 |
|----------------------------------------------------|
| 3 AB 1 A 1 0 160.7 |
| 3 AB 2 B 0 1 156.8 |
+----------------------------------------------------+
.
. * participants per sequence (one row per participant)
. egen byte first = tag(id)
. tabulate seq if first
seq | Freq. Percent Cum.
------------+-----------------------------------
AB | 12 50.00 50.00
BA | 12 50.00 100.00
------------+-----------------------------------
Total | 24 100.00
.
. * cell means by sequence and period
. table seq period, statistic(mean sbp) nformat(%9.2f)
-----------------------------------
| period
| 1 2 Total
--------+--------------------------
seq |
AB | 145.72 149.60 147.66
BA | 150.01 137.86 143.93
Total | 147.86 143.73 145.80
-----------------------------------
.
. * mixed model: treatment and period as fixed effects, a random intercept per participant
. * (REML, Kenward-Roger degrees of freedom). Fitted quietly; the lines below print its treatment
. * and period rows and the two standard deviations. R prints the same model in full.
. quietly mixed sbp i.trt i.period || id:, reml dfmethod(kroger)
. matrix T = r(table)
. scalar mm_d = T[rownumb(T, "b"), colnumb(T, "sbp:1.trt")]
. scalar mm_se = T[rownumb(T, "se"), colnumb(T, "sbp:1.trt")]
. scalar mm_lo = T[rownumb(T, "ll"), colnumb(T, "sbp:1.trt")]
. scalar mm_hi = T[rownumb(T, "ul"), colnumb(T, "sbp:1.trt")]
. scalar mm_df = T[rownumb(T, "df"), colnumb(T, "sbp:1.trt")]
. scalar mm_p = T[rownumb(T, "pvalue"), colnumb(T, "sbp:1.trt")]
. scalar mm_pi = T[rownumb(T, "b"), colnumb(T, "sbp:2.period")]
. scalar mm_pilo = T[rownumb(T, "ll"), colnumb(T, "sbp:2.period")]
. scalar mm_pihi = T[rownumb(T, "ul"), colnumb(T, "sbp:2.period")]
. scalar mm_sdid = exp(_b[lns1_1_1:_cons])
. scalar mm_sdres = exp(_b[lnsig_e:_cons])
. display "mixed model, treatment (A minus B): " %8.4f mm_d " SE " %7.4f mm_se " 95% CI " %8.4f mm_lo " to " %8.4f mm_
> hi
mixed model, treatment (A minus B): -8.0167 SE 1.7938 95% CI -11.7367 to -4.2966
. display "mixed model, treatment: df " %7.4f mm_df " P " %6.4f mm_p
mixed model, treatment: df 22.0000 P 0.0002
. display "mixed model, period (2 minus 1): " %8.4f mm_pi " 95% CI " %8.4f mm_pilo " to " %8.4f mm_pihi
mixed model, period (2 minus 1): -4.1333 95% CI -7.8534 to -0.4133
. display "mixed model, SD between participants " %7.4f mm_sdid " SD within participants " %7.4f mm_sdres
mixed model, SD between participants 5.1869 SD within participants 6.2139
.
. * mixed model with an explicit carry-over term: treatment is then estimated from period 1 only
. quietly mixed sbp i.trt i.period i.carry || id:, reml dfmethod(kroger)
. matrix C = r(table)
. scalar mmc_d = C[rownumb(C, "b"), colnumb(C, "sbp:1.trt")]
. scalar mmc_lo = C[rownumb(C, "ll"), colnumb(C, "sbp:1.trt")]
. scalar mmc_hi = C[rownumb(C, "ul"), colnumb(C, "sbp:1.trt")]
. scalar mmc_l = C[rownumb(C, "b"), colnumb(C, "sbp:1.carry")]
. scalar mmc_lse = C[rownumb(C, "se"), colnumb(C, "sbp:1.carry")]
. scalar mmc_lp = C[rownumb(C, "pvalue"), colnumb(C, "sbp:1.carry")]
. display "carry-over model, treatment (A minus B): " %8.4f mmc_d " 95% CI " %8.4f mmc_lo " to " %8.4f mmc_hi
carry-over model, treatment (A minus B): -4.2917 95% CI -10.8943 to 2.3110
. display "carry-over model, carry-over term: " %8.4f mmc_l " SE " %7.4f mmc_lse " P " %6.4f mmc_lp
carry-over model, carry-over term: 7.4500 SE 5.4483 P 0.1853
.
. * one row per participant: period difference d = period 1 minus period 2, and the participant total
. keep id seq period sbp
. reshape wide sbp, i(id) j(period)
(j = 1 2)
Data Long -> Wide
-----------------------------------------------------------------------------
Number of observations 48 -> 24
Number of variables 4 -> 4
j variable (2 values) period -> (dropped)
xij variables:
sbp -> sbp1 sbp2
-----------------------------------------------------------------------------
. generate double d = sbp1 - sbp2
. generate double tot = sbp1 + sbp2
. generate byte grp = cond(seq == "AB", 1, 2)
. label define grp 1 "AB" 2 "BA"
. label values grp grp
.
. * treatment effect from period differences: Delta-hat = (mean d in AB - mean d in BA)/2, two-sample t-test
. ttest d, by(grp)
Two-sample t test with equal variances
------------------------------------------------------------------------------
Group | Obs Mean Std. err. Std. dev. [95% conf. interval]
---------+--------------------------------------------------------------------
AB | 12 -3.883336 2.973998 10.30223 -10.42906 2.66239
BA | 12 12.15 2.006485 6.950669 7.733753 16.56624
---------+--------------------------------------------------------------------
Combined | 24 4.133331 2.423217 11.87129 -.8794746 9.146137
---------+--------------------------------------------------------------------
diff | -16.03333 3.587569 -23.4735 -8.593171
------------------------------------------------------------------------------
diff = mean(AB) - mean(BA) t = -4.4691
H0: diff = 0 Degrees of freedom = 22
Ha: diff < 0 Ha: diff != 0 Ha: diff > 0
Pr(T < t) = 0.0001 Pr(|T| > |t|) = 0.0002 Pr(T > t) = 0.9999
. scalar dbar_ab = r(mu_1)
. scalar dbar_ba = r(mu_2)
. scalar df_d = r(df_t)
. scalar delta = (r(mu_1) - r(mu_2)) / 2
. scalar delta_se = r(se) / 2
. scalar delta_p = r(p)
. scalar tc = invttail(df_d, 0.025)
. display "treatment (A minus B): " %8.4f delta " SE " %7.4f delta_se " 95% CI " %8.4f delta - tc * delta_se " to " %8
> .4f delta + tc * delta_se " P " %6.4f delta_p
treatment (A minus B): -8.0167 SE 1.7938 95% CI -11.7367 to -4.2966 P 0.0002
.
. * period effect (period 2 minus period 1) from the same differences: -(mean d in AB + mean d in BA)/2
. scalar per = -(dbar_ab + dbar_ba) / 2
. display "period (2 minus 1): " %8.4f per " 95% CI " %8.4f per - tc * delta_se " to " %8.4f per + tc * delta_se
period (2 minus 1): -4.1333 95% CI -7.8534 to -0.4132
.
. * carry-over test: compare participant totals between sequences (between-participant, so low power)
. ttest tot, by(grp)
Two-sample t test with equal variances
------------------------------------------------------------------------------
Group | Obs Mean Std. err. Std. dev. [95% conf. interval]
---------+--------------------------------------------------------------------
AB | 12 295.3167 3.709261 12.84926 287.1526 303.4807
BA | 12 287.8667 3.990715 13.82424 279.0832 296.6502
---------+--------------------------------------------------------------------
Combined | 24 291.5917 2.7752 13.59565 285.8507 297.3326
---------+--------------------------------------------------------------------
diff | 7.449999 5.448341 -3.849168 18.74917
------------------------------------------------------------------------------
diff = mean(AB) - mean(BA) t = 1.3674
H0: diff = 0 Degrees of freedom = 22
Ha: diff < 0 Ha: diff != 0 Ha: diff > 0
Pr(T < t) = 0.9073 Pr(|T| > |t|) = 0.1853 Pr(T > t) = 0.0927
.
. * design values used to simulate: sigma_s = 12, sigma_e = 6, lambda = -2
. * power of the carry-over test at lambda = -2, and the carry-over it detects with 80 percent power
. * (two-sided 5 percent level); SE of lambda-hat = sqrt((4 sigma_s^2 + 2 sigma_e^2) (1/n_AB + 1/n_BA))
. scalar sig_s = 12
. scalar sig_e = 6
. scalar lam = -2
. count if grp == 1
12
. scalar n_ab = r(N)
. count if grp == 2
12
. scalar n_ba = r(N)
. scalar des_se = sqrt((4 * sig_s^2 + 2 * sig_e^2) * (1 / n_ab + 1 / n_ba))
. scalar des_df = n_ab + n_ba - 2
. scalar des_t = invttail(des_df, 0.025)
. scalar des_power = nt(des_df, lam / des_se, -des_t) + nttail(des_df, lam / des_se, des_t)
. scalar des_mde80 = (des_t + invttail(des_df, 0.20)) * des_se
. display "design SE of the carry-over estimate: " %7.4f des_se
design SE of the carry-over estimate: 10.3923
. display "power of the carry-over test at lambda = -2: " %6.4f des_power
power of the carry-over test at lambda = -2: 0.0539
. display "carry-over detected with 80% power: " %7.4f des_mde80
carry-over detected with 80% power: 30.4717
R: the same analysis with lme4 and lmerTest
# Analysis of a 2x2 AB/BA crossover (simulated data, not evidence about any real drug)
# Outcome sbp = systolic blood pressure in mmHg; drug A versus drug B; 12 + 12 participants, two periods.
# Packages: lme4 and lmerTest, plus pbkrtest for the Kenward-Roger degrees of freedom.
suppressPackageStartupMessages(library(lmerTest))
# read the simulated dataset (not published with this article)
dat <- read.csv("crossover.csv")
str(dat)
print(head(dat, 6))
# participants per sequence (one row per participant)
ids <- unique(dat[, c("id", "seq")])
print(table(ids$seq))
# cell means by sequence and period
print(round(tapply(dat$sbp, list(seq = dat$seq, period = dat$period), mean), 2))
# mixed model: treatment and period as fixed effects, a random intercept per participant (REML, Kenward-Roger df)
dat$trt <- factor(dat$trt)
dat$period <- factor(dat$period)
dat$carry <- factor(dat$carry)
mm <- lmer(sbp ~ trt + period + (1 | id), data = dat, REML = TRUE)
smm <- summary(mm, ddf = "Kenward-Roger")
print(smm)
cm <- coef(smm)
tq <- qt(0.975, cm["trt1", "df"])
# 95% CI of the treatment coefficient (A minus B)
print(round(cm["trt1", "Estimate"] + c(lo = -1, hi = 1) * tq * cm["trt1", "Std. Error"], 4))
# mixed model with an explicit carry-over term: treatment is then estimated from period 1 only
mmc <- lmer(sbp ~ trt + period + carry + (1 | id), data = dat, REML = TRUE)
smmc <- summary(mmc, ddf = "Kenward-Roger")
print(smmc)
cc <- coef(smmc)
tqc <- qt(0.975, cc["trt1", "df"])
# 95% CI of the treatment coefficient in this model
print(round(cc["trt1", "Estimate"] + c(lo = -1, hi = 1) * tqc * cc["trt1", "Std. Error"], 4))
# one row per participant: period difference d = period 1 minus period 2, and the participant total
w <- reshape(dat[, c("id", "seq", "period", "sbp")], idvar = c("id", "seq"),
timevar = "period", direction = "wide")
w$d <- w$sbp.1 - w$sbp.2
w$tot <- w$sbp.1 + w$sbp.2
w$seq <- factor(w$seq, levels = c("AB", "BA"))
# treatment effect from period differences: Delta-hat = (mean d in AB - mean d in BA)/2, two-sample t-test
td <- t.test(d ~ seq, data = w, var.equal = TRUE)
print(td)
dbar_ab <- unname(td$estimate[1])
dbar_ba <- unname(td$estimate[2])
delta <- (dbar_ab - dbar_ba) / 2
delta_se <- td$stderr / 2
tc <- qt(0.975, unname(td$parameter))
cat(sprintf("treatment (A minus B): %.4f, SE %.4f, 95%% CI %.4f to %.4f\n",
delta, delta_se, delta - tc * delta_se, delta + tc * delta_se))
# period effect (period 2 minus period 1) from the same differences: -(mean d in AB + mean d in BA)/2
per <- -(dbar_ab + dbar_ba) / 2
cat(sprintf("period (2 minus 1): %.4f, 95%% CI %.4f to %.4f\n", per, per - tc * delta_se, per + tc * delta_se))
# carry-over test: compare participant totals between sequences (between-participant, so low power)
tt <- t.test(tot ~ seq, data = w, var.equal = TRUE)
print(tt)
# design values used to simulate: sigma_s = 12, sigma_e = 6, lambda = -2
# power of the carry-over test at lambda = -2, and the carry-over it detects with 80 percent power
# (two-sided 5 percent level); SE of lambda-hat = sqrt((4 sigma_s^2 + 2 sigma_e^2) (1/n_AB + 1/n_BA))
sigma_s <- 12; sigma_e <- 6; lambda <- -2
n_ab <- sum(ids$seq == "AB"); n_ba <- sum(ids$seq == "BA")
des_se <- sqrt((4 * sigma_s^2 + 2 * sigma_e^2) * (1 / n_ab + 1 / n_ba))
des_df <- n_ab + n_ba - 2
des_t <- qt(0.975, des_df)
des_power <- pt(-des_t, des_df, lambda / des_se) + 1 - pt(des_t, des_df, lambda / des_se)
des_mde80 <- (des_t + qt(0.80, des_df)) * des_se
cat(sprintf("design SE of the carry-over estimate: %.4f\n", des_se))
cat(sprintf("power of the carry-over test at lambda = -2: %.4f\n", des_power))
cat(sprintf("carry-over detected with 80%% power: %.4f\n", des_mde80))
> suppressPackageStartupMessages(library(lmerTest))
> dat <- read.csv("crossover.csv")
> str(dat)
'data.frame': 48 obs. of 7 variables:
$ id : int 1 1 2 2 3 3 4 4 5 5 ...
$ seq : chr "BA" "BA" "AB" "AB" ...
$ period : int 1 2 1 2 1 2 1 2 1 2 ...
$ treatment: chr "B" "A" "A" "B" ...
$ trt : int 0 1 1 0 1 0 1 0 1 0 ...
$ carry : int 0 0 0 1 0 1 0 1 0 1 ...
$ sbp : num 146 134 152 150 161 ...
> print(head(dat, 6))
id seq period treatment trt carry sbp
1 1 BA 1 B 0 0 146.2
2 1 BA 2 A 1 0 134.5
3 2 AB 1 A 1 0 152.4
4 2 AB 2 B 0 1 149.6
5 3 AB 1 A 1 0 160.7
6 3 AB 2 B 0 1 156.8
> ids <- unique(dat[, c("id", "seq")])
> print(table(ids$seq))
AB BA
12 12
> print(round(tapply(dat$sbp, list(seq = dat$seq, period = dat$period),
+ mean), 2))
period
seq 1 2
AB 145.72 149.60
BA 150.01 137.86
> dat$trt <- factor(dat$trt)
> dat$period <- factor(dat$period)
> dat$carry <- factor(dat$carry)
> mm <- lmer(sbp ~ trt + period + (1 | id), data = dat,
+ REML = TRUE)
> smm <- summary(mm, ddf = "Kenward-Roger")
> print(smm)
Linear mixed model fit by REML. t-tests use Kenward-Roger's method [
lmerModLmerTest]
Formula: sbp ~ trt + period + (1 | id)
Data: dat
REML criterion at convergence: 321
Scaled residuals:
Min 1Q Median 3Q Max
-2.0487 -0.4621 -0.1869 0.6003 1.4972
Random effects:
Groups Name Variance Std.Dev.
id (Intercept) 26.90 5.187
Residual 38.61 6.214
Number of obs: 48, groups: id, 24
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 151.871 1.880 44.797 80.784 < 2e-16 ***
trt1 -8.017 1.794 22.000 -4.469 0.000192 ***
period2 -4.133 1.794 22.000 -2.304 0.031028 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) trt1
trt1 -0.477
period2 -0.477 0.000
> cm <- coef(smm)
> tq <- qt(0.975, cm["trt1", "df"])
> print(round(cm["trt1", "Estimate"] + c(lo = -1, hi = 1) *
+ tq * cm["trt1", "Std. Error"], 4))
lo hi
-11.7367 -4.2966
> mmc <- lmer(sbp ~ trt + period + carry + (1 | id),
+ data = dat, REML = TRUE)
> smmc <- summary(mmc, ddf = "Kenward-Roger")
> print(smmc)
Linear mixed model fit by REML. t-tests use Kenward-Roger's method [
lmerModLmerTest]
Formula: sbp ~ trt + period + carry + (1 | id)
Data: dat
REML criterion at convergence: 313.9
Scaled residuals:
Min 1Q Median 3Q Max
-2.1949 -0.5026 -0.0767 0.4785 1.5487
Random effects:
Groups Name Variance Std.Dev.
id (Intercept) 25.22 5.022
Residual 38.61 6.214
Number of obs: 48, groups: id, 24
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 150.008 2.306 38.059 65.041 <2e-16 ***
trt1 -4.292 3.262 38.059 -1.316 0.1961
period2 -7.858 3.262 38.059 -2.409 0.0209 *
carry1 7.450 5.448 22.000 1.367 0.1853
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) trt1 perid2
trt1 -0.707
period2 0.279 -0.698
carry1 -0.591 0.835 -0.835
> cc <- coef(smmc)
> tqc <- qt(0.975, cc["trt1", "df"])
> print(round(cc["trt1", "Estimate"] + c(lo = -1, hi = 1) *
+ tqc * cc["trt1", "Std. Error"], 4))
lo hi
-10.8943 2.3110
> w <- reshape(dat[, c("id", "seq", "period", "sbp")],
+ idvar = c("id", "seq"), timevar = "period", direction = "wide")
> w$d <- w$sbp.1 - w$sbp.2
> w$tot <- w$sbp.1 + w$sbp.2
> w$seq <- factor(w$seq, levels = c("AB", "BA"))
> td <- t.test(d ~ seq, data = w, var.equal = TRUE)
> print(td)
Two Sample t-test
data: d by seq
t = -4.4691, df = 22, p-value = 0.0001918
alternative hypothesis: true difference in means between group AB and group BA is not equal to 0
95 percent confidence interval:
-23.473498 -8.593169
sample estimates:
mean in group AB mean in group BA
-3.883333 12.150000
> dbar_ab <- unname(td$estimate[1])
> dbar_ba <- unname(td$estimate[2])
> delta <- (dbar_ab - dbar_ba)/2
> delta_se <- td$stderr/2
> tc <- qt(0.975, unname(td$parameter))
> cat(sprintf("treatment (A minus B): %.4f, SE %.4f, 95%% CI %.4f to %.4f\n",
+ delta, delta_se, delta - tc * delta_se, delta + tc * delta_se))
treatment (A minus B): -8.0167, SE 1.7938, 95% CI -11.7367 to -4.2966
> per <- -(dbar_ab + dbar_ba)/2
> cat(sprintf("period (2 minus 1): %.4f, 95%% CI %.4f to %.4f\n",
+ per, per - tc * delta_se, per + tc * delta_se))
period (2 minus 1): -4.1333, 95% CI -7.8534 to -0.4133
> tt <- t.test(tot ~ seq, data = w, var.equal = TRUE)
> print(tt)
Two Sample t-test
data: tot by seq
t = 1.3674, df = 22, p-value = 0.1853
alternative hypothesis: true difference in means between group AB and group BA is not equal to 0
95 percent confidence interval:
-3.849168 18.749168
sample estimates:
mean in group AB mean in group BA
295.3167 287.8667
> sigma_s <- 12
> sigma_e <- 6
> lambda <- -2
> n_ab <- sum(ids$seq == "AB")
> n_ba <- sum(ids$seq == "BA")
> des_se <- sqrt((4 * sigma_s^2 + 2 * sigma_e^2) * (1/n_ab +
+ 1/n_ba))
> des_df <- n_ab + n_ba - 2
> des_t <- qt(0.975, des_df)
> des_power <- pt(-des_t, des_df, lambda/des_se) + 1 -
+ pt(des_t, des_df, lambda/des_se)
> des_mde80 <- (des_t + qt(0.8, des_df)) * des_se
> cat(sprintf("design SE of the carry-over estimate: %.4f\n",
+ des_se))
design SE of the carry-over estimate: 10.3923
> cat(sprintf("power of the carry-over test at lambda = -2: %.4f\n",
+ des_power))
power of the carry-over test at lambda = -2: 0.0539
> cat(sprintf("carry-over detected with 80%% power: %.4f\n",
+ des_mde80))
carry-over detected with 80% power: 30.4717
Pitfalls
-
Crediting a period-2 change to the period-2 drug
The change from period 1 to period 2 inside one sequence is read as the effect of the drug taken in period 2, although part of it may belong to period 2 itself.
Fix: Compare the sequences, not the periods within one sequence. Half the difference between the two mean period differences removes any shift shared by period 2.
-
Treating the sequence effect as a separate effect that can be estimated
In a 2x2 crossover, a difference between the AB and BA groups cannot be split into sequence, carry-over and treatment-by-period interaction.
Fix: Report the comparison of totals as one aliased quantity, and prevent carry-over by design.
-
"A non-significant carry-over test shows there is no carry-over."
In the simulated trial the test gives P = 0.19, although a carry-over of -2 mmHg was built into the data.
Fix: The carry-over test compares participant totals between sequences, a between-participant comparison with little power. A non-significant result does not show that carry-over is absent; carry-over is prevented by an adequate washout, not removed by a preliminary test.
-
Testing carry-over first, then falling back to period 1
The two-stage procedure inflates the type I error and biases the estimate of the treatment effect [5].
Fix: Prespecify the period-difference analysis, or the equivalent mixed model, and plan an adequate washout.
What to do in your own analysis
- Plan carry-over out at the design stage, with a stable condition and a washout that the drug's known duration of action suggests is long enough.
- Consider prespecifying the primary analysis as the period-difference estimate or, equivalently with complete data, a mixed model with treatment, period and a random intercept per participant.
- Report the period estimate beside the treatment estimate, so readers can see how much time alone moved the outcome.
- If a carry-over estimate is reported, present it as a low-powered between-participant comparison that does not choose the analysis.
- Check the report against the CONSORT extension for randomised crossover trials, the CONSORT reporting guideline adapted to this design [6]. Among other items, it asks for the rationale for a crossover design and for statistical methods that take the paired nature of the data into account.
Glossary
- crossover trial (การศึกษาแบบไขว้)
- A trial in which each participant receives every treatment under study, one per period, in a randomised order.
- sequence (ลำดับ)
- The order in which a participant receives the treatments; in a 2x2 crossover, AB or BA.
- period effect (ผลของช่วงเวลา)
- A shift in the outcome shared by all participants in one period compared with another, whichever treatment they take.
- carry-over effect (ผลตกค้าง)
- An effect of an earlier period's treatment that persists into a later period.
- washout period (ช่วงล้างยา)
- A treatment-free interval between periods that lets the first treatment's effect fade.
- aliasing (การปนกันของผล)
- Two effects are aliased when the design gives them the same pattern in the data, so no analysis can separate them.
- within-participant comparison (การเปรียบเทียบภายในผู้เข้าร่วม)
- A comparison of readings taken on the same person, so stable differences between people cancel.
- power (อำนาจการทดสอบ)
- The probability that a test declares an effect when an effect of a stated size truly exists.
- CONSORT extension (ส่วนขยายของ CONSORT)
- A version of the CONSORT reporting guideline for randomised trials adapted to one design; the crossover version appeared in 2019.
References
- Senn S. Cross-over trials in clinical research. 2nd ed. Wiley; 2002. https://doi.org/10.1002/0470854596
- Jones B, Kenward MG. Design and analysis of cross-over trials. 3rd ed. Chapman and Hall/CRC; 2014. https://doi.org/10.1201/b17537
- Hills M, Armitage P. The two-period cross-over clinical trial. Br J Clin Pharmacol. 1979;8(1):7-20. https://doi.org/10.1111/j.1365-2125.1979.tb05903.x
- Grizzle JE. The two-period change-over design and its use in clinical trials. Biometrics. 1965;21(2):467-80. https://doi.org/10.2307/2528104
- Freeman PR. The performance of the two-stage analysis of two-treatment, two-period crossover trials. Stat Med. 1989;8(12):1421-32. https://doi.org/10.1002/sim.4780081202
- Dwan K, Li T, Altman DG, Elbourne D. CONSORT 2010 statement: extension to randomised crossover trials. BMJ. 2019;366:l4378. https://doi.org/10.1136/bmj.l4378
Key takeaways
- In a 2x2 crossover, half the difference between the two sequences' mean period differences estimates the treatment effect, and any period effect cancels.
- With complete data, a random-intercept mixed model with treatment and period estimates half the difference in mean period differences, with the P value and degrees of freedom of the pooled-variance t-test.
- Carry-over is aliased with sequence: it shifts the treatment estimate by minus half its size and can be estimated only from comparisons between participants.
- The carry-over test has little power, so a non-significant result does not show that carry-over is absent.
- Prevent carry-over with an adequate washout, and do not let a preliminary carry-over test choose the analysis.
Related in the wiki: [[treatment-by-time-interaction-longitudinal-trials]] [[crossover-trials-design-analysis-guide]] [[mixed-model-series-guide]]