クロスオーバー試験は、生物学的同等性試験や薬物動態試験で利用されています。ウォッシュアウト期間の設定がとても大切で、治療効果が無くなるのを十分に待つことになります。このデザインはリハビリテーション研究でも利用されていますが、持ち越し効果や時期効果が十分に吟味されず、間違った使い方をしているものも見られます。これらの効果はP値に直接関連してきますので、十分に注意しましょう。
ここでは、折笠(2016)のなかにあるサンプルデータを使用して、同じ解析結果をRで算出してみます。FEVの値をpointとして、有意水準は5％で検定します。
クロスオーバー試験の研究デザイン
2群2期のクロスオーバーデザイン、2×2 デザイン、AB/BA デザイン と呼ばれています
ここでは $ Y_{ij} $ という標記を使います(2×２の行列と考えてください)
$\blacklozenge \quad i $ は順序群（行）を示す ➡ 1＝Sequence1 (AB群)、2＝Sequence2 (BA群)
$\blacklozenge \quad j $ は時期（列）を示す ➡ 1＝I期、2＝II期
$\blacklozenge \quad Y_{ijk} $ は実測値を示す
データの準備と要約
使用するパッケージ（パッケージのインストール）
library(lmerTest)
library(exactRankTests)
データセット crossover をdatに格納します（ファイルの読み込み）
dat <- read.csv("crossover.csv", header=T, fileEncoding = "UTF-8")
Y11, Y21, Y12, Y22のデータの要約
#IDを名義尺度に変換
dat$ID <- factor(dat$ID)
参考文献に合わせて基準を変更します
基準を変更する場合にはfactor変数にしなければなりません
Rは数字であれば若い数、アルファベットはAに近い文字が基準となります。差の検定を行うと、II－IやB－Aの解を出力します。折笠（2016）ではI－IIやA－Bの解を出力しています。factor変数に変換して、levels関数で調べたら分かります。
もちろんChat GPTでも教えてくれますよ
#factor変数に変更
dat$period <- factor(dat$period)
dat$treatment <- factor(dat$treatment)
#基準を調べてみます（最初に書いてある文字が基準です）
levels(dat$period)
levels(dat$treatment)
基準を変更します
dat$period <- relevel(dat$period, "II")
dat$treatment <- relevel(dat$treatment, "B")
#確認
levels(dat$period)
levels(dat$treatment)
#I期
#sequence1の治療A
Y11 <- dat[dat$sequence=="s1" & dat$treatment=="A", "point"]
#sequence2の治療B
Y21 <- dat[dat$sequence=="s2" & dat$treatment=="B", "point"]
#II期
#sequence1の治療B
Y12 <- dat[dat$sequence=="s1" & dat$treatment=="B", "point"]
#sequence2の治療A
Y22 <- dat[dat$sequence=="s2" & dat$treatment=="A", "point"]
#データの要約
summary(Y11)
summary(Y21)
summary(Y12)
summary(Y22)
箱ひげ図
name <- c("Y11", "Y21", "Y12", "Y22")
boxplot(Y11, Y21, Y12, Y22, xaxt="n")
axis(side=1, at=1:4, labels=name)
参考文献（折笠 2016）の図3の描き方
fig <- function(){
plot(
c(1, 1.2, 1.8, 2),
c(mean(Y11), mean(Y21), mean(Y12), mean(Y22)),
xlab = "", ylab = "",
ylim = c(0,5), xlim = c(0.5, 2.5),
xaxt="n")
name <- c("I期", "II期")
axis(side=1, at=c(1.1, 1.9), labels=name)
}
fig()
#標準誤差の関数を作成
stm <- function(x) sd(x)/sqrt(length(x))
arrows(
c(1, 1.2, 1.8, 2),
c(
mean(Y11) - stm(Y11),
mean(Y21) - stm(Y21),
mean(Y12) - stm(Y12),
mean(Y22) - stm(Y22)
),
c(1, 1.2, 1.8, 2),
c(
mean(Y11) + stm(Y11),
mean(Y21) + stm(Y21),
mean(Y12) + stm(Y12),
mean(Y22) + stm(Y22)
),
code=3, angle=90
)
クロスオーバー試験で想定している結果
対象者を各順序群に無作為に割り付けることで、両群の治療Aと治療Bの差が同じくらいになることを想定しています
効果の検証
持ち越し効果（carry-over effect）がある場合
washout期間は十分に確保できていたか？
＊後発医薬品の生物学的同等性試験ガイドラインでは、washout期間 (休薬期間) が設定されています
持ち越し効果の検定
各sequenceにおける治療A＋治療Bの平均値の差の検定（対応のないt検定）により検証
検出力が低いので有意水準は10％が良い（折笠, 2016）
検定の結果 p<0.1 になった場合は、クロスオーバー試験は適用できません
$\blacklozenge$ サンプルの持ち越し効果の検定
$Y_{11k}+Y_{12k}$ と $Y_{21k}+Y_{22k}$ の平均の差の検定
参考文献（折笠 2016）の図4の結果
#ｔ検定
t1 <- t.test(c(Y11 + Y12), c(Y21 + Y22) , var.equal=T)
#ｔ検定の結果
print(t1)
#AB群平均－BA群平均
mean(Y11 + Y12) - mean(Y21 + Y22)
#AB群平均－BA群平均の標準誤差
t1$stderr
p=0.1253（>0.1）のため持ち越し効果はないと判定します
治療効果（treatment effect）の検証 ＜ｔ検定と混合効果モデル＞
対応のあるt検定
治療Aと治療Bの比較
参考文献（折笠 2016）の図8の結果
#対応のあるｔ検定（治療A－治療B）
t2 <- t.test(c(Y11, Y22) - c(Y12, Y21) )
#ｔ検定の結果
print(t2)
#AB群平均－BA群平均の標準誤差
t2$stderr
p=0.042なので治療Bが有意に高いことになります
混合効果モデルでも対応のあるt検定と同じ解を得ることができます
t3 <- lmerTest::lmer(
point ~ treatment + (1|ID),
data = dat
)
summary(t3)
A－B＝－0.2647, t値=－2.21, p値＝0.0421
ウィルコクソンの符号順位検定では以下のような結果となります (パッケージexactRankTestsを使用)
exactRankTests::wilcox.exact(
c(Y11, Y22) , c(Y12, Y21),
paired=T, conf.int=T
)
検定統計量の結果が参考文献（折笠 2016）の図8の結果と異なります。詳細が分かれば追加で記載します。
混合効果モデル（対象者ID＝変量効果）
参考文献（折笠 2016）の図5の結果
period effectを考慮した治療Aと治療Bの差の検定
対応のある検定になるので対象者（ID）を変量効果とした混合効果モデルを適用
fit1 <- lmerTest::lmer(
point ~ treatment + period + (1|ID),
data = dat
)
summary(fit1)
period効果を考慮した場合の治療A－治療B平均値の推定値は－0.2565
p値＝0.0472なので治療Bが有意に高いことが分かります
時期効果（period effect）の検証
参考文献（折笠 2016）の図6の結果
period と treatment の位置を変えてますが、上述した混合効果モデルと同じです
fit2 <- lmer(
point ~ period + treatment + (1|ID),
data = dat
)
summary(fit2)
treatmentで調整したI期—II期の平均値の推定値は0.139
p=0.2595なので有意差は認められません
反復測定分散分析：repeated measures ANOVA （混合効果モデル）
ID（被験者）を変量効果とした混合効果モデル
参考文献（折笠 2016）の図10の結果
$\blacklozenge$ 固定効果の検定
fit3 <- lmer(
point ~ sequence + treatment + period + (1|ID),
data = dat
)
anova(fit3)
平方和は分散分析表から求めることができます
$\blacklozenge$ 分散分析表 (Between Subject) の平方和
Between SubjectですのでIDを説明変数に投入した線形回帰モデル
fit4 <- lm(
point ~ sequence + ID,
data = dat
)
anova(fit4)
$\blacklozenge$ 分散分析表 (Within-subject) の平方和
Within-subjectですのでIDを変量効果とした混合効果モデル
fit5 <- lmer(
point ~ treatment + period + (1|ID),
data = dat
)
anova(fit5)
事後検定
参考文献（折笠 2016）の図11の結果
分散分析を適用します
$\blacklozenge$ 効果の検定
fit6 <- lm(
point ~ period + treatment + ID,
data = dat
)
anova(fit6)
periodの結果が参考文献（折笠 2016）の図11の結果と異なります。原因が分かり次第修正します。
最小二乗平均
lmerTestのls_means関数を使用
ls_means(fit5)
$\blacklozenge$ I期、II期、治療A、治療Bの平均
mean(dat[dat$period=="I", "point"])
mean(dat[dat$period=="II", "point"])
mean(dat[dat$treatment=="A", "point"])
mean(dat[dat$treatment=="B", "point"])
最小二乗平均差
lmerTestのdifflsmeans関数
#参考文献に合わせます
dat$period <- relevel(dat$period, "I")
dat$treatment <- relevel(dat$treatment, "A")
#確認
levels(dat$period)
levels(dat$treatment)
fit5_2 <- lmer(
point ~ treatment + period + (1|ID),
data = dat
)
difflsmeans(fit5_2)
求め方・・・
$\blacklozenge$ 各群の最小二乗法推定量
繰り返しありと仮定して二元配置分散分析から求めます
fit7 <- lm(
point ~ period + treatment,
data=dat
)
summary(fit7)$coef
Y11（IのA）：1.9611029 + 0.1390278 – 0.2565278 = 1.843603
Y12（IIのB）：1.9611029
Y21（IのB）：1.9611029 + 0.1390278 = 2.100131
Y22（IIのA）：1.9611029 – 0.2565278 = 1.704575
I期：(1.8436029+2.100131)/2=1.971867
II期：(1.961103+1.704575)/2=1.832839
治療A：(1.8436029+1.704575)/2=1.774089
治療B：(1.961103+2.100131)/2=2.030617
$\blacklozenge$ 差の平均の最小二乗法推定量
summary(fit6)$coef[c("periodI", "treatmentA"), ]
I期－II期：0.1390278
治療A－治療B：－0.2565278
$\blacklozenge$ 差の平均の最小二乗法推定量の95%信頼区間
confint ( fit6 , level = 0.95 )[c("periodI", "treatmentA"), ]
参考文献
折笠秀樹. クロスオーバー試験の計画および解析. 薬理と治療, 2016, 44.9: 1261-1276.
角間辰之, 服部聡; 臨床試験のデザインと解析: 薬剤開発のためのバイオ統計. 近代科学社. 2012
コメント欄 『間違い』や『分かりにくい部分』などのご意見もお寄せください