Document

SASによる二項比率の差の非劣性を示す
ベイズ指標の算出方法について
川﨑洋平 1), 2)
榊原伊織 2)
宮岡悦良 3)
1) 国立国際医療研究センター 2) タクミインフォメーションテクノロジー 3) 東京理科大学
On calculation method of a Bayesian noninferiority index for two binomial proportions
using SAS
Yohei Kawasaki 1), 2) Iori Sakakibara 2) Etsuo Miyaoka 3)
1) National Center for Global Health and Medicine
2) Takumi Information Technology Inc.
3) Tokyo University of Science
E-mail: [email protected] URL:https://sites.google.com/site/biostatncgmjp/welcome
要旨：
Kawasaki and Miyaoka (2013)は二項比率の差の非劣性を示
す指標
τ = P( π1 > π2 – Δ0 | X1, X2)
を提案した．
本研究では，この指標を紹介するとともに，FCMP Procedure
を用いた，SASプログラムを提供する．
キーワード：
Non-inferiority Test, Bayes Inference, Appell Function,
FCMP Procedure
2
本日の発表構成
□ 本日の発表構成
1. 母比率の差の非劣性を表す指標 τ の紹介
τ = P( π1 > π2 － Δ0 | x1, x2 )
（π1: 治験薬の母比率，π2：対照薬の母比率，Δ0：非劣性マージン 0 < Δ0< 1 ）
⇒ 治験薬の母比率が対照薬の母比率よりΔ0以上劣らない確率
2. τ に関する計算式の紹介
2.1 近似的に求める方法及びプログラムの紹介（近似確率）
2.2 正確に求める方法及びプログラムの紹介（正確確率）
3. 実際の臨床試験への適用
4. Macroの紹介とまとめ
3
臨床試験の例
審査報告書より
□ ゴナールエフ vs フェルティノームのPhase-3 試験
対象疾患：無排卵周期症
治験薬：ゴナールエフ
対照薬：フェルティノーム
□ 試験結果
薬剤
治験薬
対照薬
排卵
あり
102
109
なし
27
23
計
129
132
□ 検証したい仮説
治験薬の改善率は対照薬の改善率より10%以上劣っていない．
□ 解析方法
Test : Dunnett and Gent (1977), Farrington and Manning (1990)
C.I. : Newcombe (1998), Agresti and Min (2001)
⇒ このようにFrequentistの枠組みで解析することが一般的である．
4
母比率の差の非劣性に関して
Frequentist and Bayesian Framework
□ Frequentist
1. 差が無いという仮定で、観測値以上に極端なデータが得られる確率を算出
2. 片側有意水準 0.025と比較し、結論を出す。
□ Bayesian
1.治験薬の母比率が対照薬よりΔ0以上劣らない確率を直接的に算出する
例)
  P(1   2  0 | x1 , x2 )  0.95
（π1: 治験薬の母比率，π2：対照薬の母比率, Δ0：非劣性マージン 0 < Δ0< 1 ）
□ Key Point
ベイズ流のFrameworkは直感的に理解がしやすい場合もある。
5
Objective
□ 本研究の目的
1. 母比率の差の非劣性を直接的に表す指標を紹介する．
  P(1   2  0 | x1 , x2 )
（π1: 治験薬の母比率, π2：対照薬の母比率, Δ0：非劣性マージン 0 < Δ0< 1 ）
2. Kawasaki and Miyaoka (2013)では提案することができなかった、SASでの計算プログラ
ムを提供する。
3. 作成したSASプログラムを用いて実際の臨床試験にτを適用してその有用性を示す．
6
Notation 1
□ 確率変数 X1 ， X2
X1～Bin(n1 , 1 ) , X 2～Bin(n2 ,  2 )
（但し、Bin ( ni , πi ) はパラメータ ni , πi の二項分布を表す。 i =1, 2）
□ 事前分布
共役事前分布：ベータ分布
1～Beta(1 , 1 ) ,  2～Beta( 2 ,  2 )
（但し、Beta ( αi , βi )はパラメータ αi , βi のベータ分布を表す。）
7
Notation 2
□ πi の事後分布
 ni  xi
1
  i (1   i ) ni  xi
 i i 1 (1   i )  i 1
x
B( i ,  i )
g i ( i | xi )   i 
1 n 
1
 i  ixi (1   i ) ni  xi
 i i 1 (1   i )  i 1d i
0 xi
B( i ,  i )
 


1
 ixi 1 1 (1   i ) ni  xi   i 1
B( xi   i , ni  xi   i )
（但し、B ( a , b ) はベータ関数を表し、i =1, 2とする）
□ Posterior Notation
1, post～Beta(a1 , b1 ) ,  2, post～Beta(a2 , b2 )
（但し、 ai  xi   i , bi  ni  xi  i とする）
8
1. 近似確率
□ 期待値・分散
E ( 1, post   2, post )  1, post   2, post
V ( 1, post   2, post ) 
1, post (1  1, post )
a1  b1  1

 2, post (1   2, post )
a2  b2  1
（但し、 i , post  ai /(ai  bi ) とする。）
□ 近似確率の算出式
  P( 1   2   0 | x1 , x2 )
  0  ( 1, post   2, post )
 1  (
1, post (1  1, post )
a1  b1  1

 2, post (1   2, post )
)
a2  b2  1
（但し、Φ(.)は標準正規分布の累積分布関数をあらわす）
9
近似確率のSASプログラム
□ プログラム
CDF関数を使用して容易に算出できる
data WORK.OUT ;
post1 = &a1 / ( &a1 + &b1 ) ;
post2 = &a2 / ( &a2 + &b2 ) ;
deno1 = ( post1 * ( 1 - post1 ) ) / ( &a1 + &b1 + 1 ) ;
deno2 = ( post2 * ( 1 - post2 ) ) / ( &a2 + &b2 + 1 ) ;
random = ( (-1) * &delta - ( post1 - post2 ) ) / sqrt( deno1 + deno2 ) ;
Tau_a = 1 - round( cdf( 'NORMAL' , random , 0 , 1 ) , 0.001 ) ;
keep Tau_a ;
run ;



  0  ( 1, post   2, post )
 a  1  
 1, post (1  1, post )  2, post (1   2, post )


a

b

1
a2  b2  1
1
1




,



10
2. 正確確率
□ Exact Probability
1
  P( 1   2   0 | x1 , x2 )  
 0
f ( )d
（但し、 f (δ) は事後分布の差（ δ = π1,post－ π2,post ）の確率密度関数を表す）
□ 事後分布の差の確率密度関数
Kawasaki and Miyaoka (2010)
 B(a2 , b1 )(1   ) a2 b1 1
F3 (a2 , b1 ,1  b2 ,1  a1 ; a2  b1 ;1   ,1   ) if

 B(a1 , b1 ) B(a2a ,bb2 )1
 B(a1 , b2 )(1   ) 1 2
f ( )  
F3 (a1 , b2 ,1  b1 ,1  a2 ; a1  b2 ;1   ;1   ) if
B
(
a
,
b
)
B
(
a
,
b
)
1 1
2 2

 B(a1  a2  1, b1  b2  1) if
 B(a , b ) B(a , b )
1 1
2 2

(但し、 F3 (k1 , k 2 , l1 , l2 ; h; u1 , u2 ) 

0   1
1    0
  0　
(k1 ) i (k 2 ) j (l1 ) i (l2 ) j u1i u2j
, max{|u1 |, | u2 | 1} )
( h) i  j
i! j!
j 0

i 0
11
正確確率のSASプログラム 1/6
□ ホッポハンマ関数
正確確率にはAPPELL関数が含まれており、その中にはホッポハンマ関数が含ま
れている。
このホッポハンマ関数を作成するところが苦労した点である。
*ポッホハンマー関数*;
* CREATE FUNCTIONS *;
options cmplib = WORK.FUNC mprint ;
proc fcmp outlib = WORK.FUNC.SAMPLE ;
* CREATE POCHHAMMER FUNCTION *;
function POCHHAMMER( val , ind ) ;
if ind = 0 then RES = 1 ;
else RES = val ;
ホッポハンマ関数をFCMPプロシジャ
do j = 1 to ind - 1 ;
を用いて, POCHHAMMER関数を作成した.
RES = RES * ( val + j ) ;
end ;
return( RES ) ;
endsub ;
12
FCMPプロシジャについて
□ FCMPプロシジャ
自身で作成した関数や数式などを関数化できる．
options cmplib = ライブラリ.SASデータセット名 ;
proc fcmp outlib =ライブラリ.SASデータセット名.エントリ名 ;
関数の情報を保存した
データセットを事前に指定
function 関数名( 引数 ) ;
if ind = 0 then RES = 1 ;
else
RES = val ;
do j = 1 to ind - 1 ;
RES = RES * ( val + j ) ;
end ;
return( RES ) ;
endsub ;
run ;
関数についての具体的な
処理の内容を記載する
戻り値を指定
関数の定義の終了
13
正確確率のSASプログラム 3/6
*ポッホハンマー関数*;
* CREATE FUNCTIONS *;
options cmplib = WORK.FUNC mprint ;
proc fcmp outlib = WORK.FUNC.SAMPLE ;
* CREATE POCHHAMMER FUNCTION *;
function POCHHAMMER( val , ind ) ;
if ind = 0 then RES = 1 ;
else RES = val ;
do j = 1 to ind - 1 ;
RES = RES * ( val + j ) ;
end ;
return( RES ) ;
endsub ;
ポッホハンマー関数：（ l1)i
l1 < 0の場合ガンマ関数が使用できないため
左のようなプログラムにて算出する
F3 (k1 , k2 , l1 , l2 ; h; u1 , u2 ) 


(k1 )i (k 2 ) j (l1 )i (l2 ) j u1i u2j
,
(
h
)
i
!
j
!
i j
j 0

i 0
14
正確確率のSASプログラム 4/6
*FDELTA関数の作成*;
同一プロシジャ内で作成した
関数を使用することができる
function FDELTA( i , j , delta , a1 , a2 , b1 , b2 ) ;
do m = 0 to i ;
do n = 0 to j ;
_L1 = POCHHAMMER( 1 - b1 , m ) ;
_L2 = POCHHAMMER( 1 - a2 , n ) ;
_K1 = gamma( a1 + m ) / gamma( a1 ) ; _K2 = gamma( b2 + n ) / gamma( b2 ) ;
_D = ( 1 + delta ) ** m ; _E = ( 1 + delta ) ** n ; _U1 = a1 + b2 + m + n ; _U2 = a1 + b2 ;
_K = ( beta( a1 , b2 ) * ( ( 1 + delta ) ** ( a1 + b2 - 1 ) ) ) /
( beta( a1 , b1 ) * beta( a2 , b2 ) ) ;
if _K1 = 0 or _K2 = 0 or _L1 = 0 or _L2 = 0 or _D = 0 or _E = 0 then KAI = 0 ;
else KAI = ( _K1 * _K2 * _L1 * _L2 * _D * _E ) /
( gamma( _U1 ) / gamma( _U2 ) * fact( n ) * fact( m ) ) ;
ADD + KAI ;
PDF = ADD * _K ;
end ;
 
(k1 )i (k 2 ) j (l1 )i (l2 ) j u1i u2j
end ;
F3 (k1 , k2 , l1 , l2 ; h; u1 , u2 ) 
,
(
h
)
i
!
j
!
return( PDF ) ;
i j
i 0 j 0
endsub ;
run ;

15
正確確率のSASプログラム 5/6
□ 確率計算
先ほどのFCMPプロシジャを用いて、作成したFDELTA関数を用いて、確率密度関数を算
出する。
data WORK._FDEL ;
do delta = -0.999 to -&delta by 0.001 ;
PDF = FDELTA( &Loop , &Loop , delta , &a1 , &a2 , &b1 , &b2 ) ; output ;
end ;
run ;
非劣性マージンまでの確率密度
関数を，FDELTA関数により算出
data WORK.OUT2 ;
set WORK._FDEL end = EOF ;
LAGPDF = LAG( PDF ) ; LAGDEL = LAG( delta ) ;
if _N_ ^= 1 then DAI = sum( LAGPDF , PDF ) * ( delta - LAGDEL ) / 2 ;
RDAI + DAI ;
Tau_e = round( 1 - RDAI , 0.001 ) ;
if EOF ;
台形法にてPDFまでの面積を算出する
run ;
最後にPDFまでの面積を1から
引くことで τ を算出する
16
正確確率のSASプログラム 6/6
□ マクロプログラムの実行結果（例）
サンプルサイズ10でハイパーパラメータをそれぞれ2, 3. 4, 5として，
非劣性マージンを0.1，繰り返し数を100回とした時
%noninf( n1 = 10 , n2 = 10 , x1 = 1 , x2 = 1 ,
alpha1 = 2 , alpha2 = 3 , beta1 = 4 , beta2 = 5 , delta = 0.1 , loop = 100 )
n1
: Sample Size 1
n2
: Sample Size 2
x1
: Outcome 1
x2
: Outcome 2
Alpha1 : Hyperparameter_α1
Alpha2 : Hyperparameter_α2
beta1 : Hyperparameter_β1
beta2 : Hyperparameter_β2
delta : Non-Inferiority Margin
loop : The number of simulations in the group
17
正確確率のSASプログラム 2/6
□ 出力結果
Item
Sample Size 1
Sample Size 2
Outcome 1
Outcome 2
Hyperparameter_α1
Hyperparameter_α2
Hyperparameter_β1
Hyperparameter_β2
Non-Inferiority Margin
Approximate Probability
Exact Probability
Value
10
10
1
1
2
3
4
5
0.10
0.686
0.690
18
Example 1
□ TJN318 vs ビフォナゾールの臨床試験
対象疾患：皮膚真菌症（カンジダ性間擦疹）
治験薬： TJN318クリーム
対照薬：ビフォナゾールクリーム (BFZ)
□ 検証仮説
H 0 : 1   2  0.1, H1 : 1   2  0.1
治験薬の改善率は対照薬の改善率より10%以上は劣っていない．
□ 解析結果
薬剤
判定
改善
非改善
計
TJN-318
39
4
43
BFZ
37
4
41
p- value
τ
Approximate
Exact
0.068
0.932
0.943
□ 結論
・ p-value = 0.068から帰無仮説は棄却できない．
・ 94.3%の確率で10%以上は劣っていないといえる．
19
Example 2
審査報告書より
□ ゴナールエフ vs フェルティノームのPhase-3 試験
対象疾患：無排卵周期症
治験薬：ゴナールエフ
対照薬：フェルティノーム
□ 検証仮説
H 0 : 1   2  0.1, H1 : 1   2  0.1 治験薬の排卵率は対照薬に10%以上劣っていない．
□ 解析結果
薬剤
治験薬
対照薬
排卵
あり
102
109
なし
27
23
計
p- value
129
132
0.093
τ
(正確確率)
0.911
□ 結論
・ p-value = 0.093から帰無仮説は棄却できない．
・ 91.1%の確率で10%以上は劣っていないといえる．
□ Key Point: Phase-3
Empirical Bayes Methodを用いて，過去の臨床試験の結果からハイパーパラメータを推定
し，再解析を行う事が出来る．
20
Empirical Bayes Method and Past Clinical Results
□ Empirical Bayes
ハイパーパラメータを以下の 2つの試験結果から推定
Marginal Likelihood
m( i ,  i | xi ) 
 n  B( xi   i , ni  xi   i )
L( i | xi ) g i ( i |  i ,  i )d i   i 
0
B( i ,  i )
 xi 
1

 m( , 
(ˆ i , ˆi )  arg max{
 i ,i
2
i
i
| xij )}
j 1
(但し、i =1,2とする。)
□ Past-Result
薬剤
試験内容
治験薬
Phase-2 国内
海外 A
対照薬
海外 A
海外 B
排卵
計
あり
なし
102
60
27
23
129
83
67
50
31
24
98
74
審査報告書より
21
□ Hyper-Parameterの推定結果
過去の2試験から推定されたHyper-Parameter
薬剤
治験薬
対照薬
Hyper Parameter
αi
βi
74.45
14.61
80.54
22.40
Prior distribution
Mean
S.D.
0.836
0.039
0.782
0.041
□ 解析結果
過去の試験結果の情報を持ったハイパーパラメータを利用して再解析
τ
排卵
薬剤
治験薬
対照薬
あり
なし
102
109
27
23
計
α i =β i =1
(無情報事前分布)
推定された
ハイパーパラメータ
129
132
0.911
0.997
□ Key Points
・治験薬の排卵率は対照薬の排卵率に10%以上は劣っていない確率が91.1% → 99.7%
となった．
・このように τ + Empirical Bayesで、過去の試験結果の情報を加味して解析が実施できる.
22
/*
+--------------------------------------------------------------------------------------------------------------------+
| DESCRIPTION: Takumi Information Technology
| VERSION
: SAS 9.2
| FILENAME : Noninf.sas
| PURPOSE
: On calculation method of a Bayesian non-inferiority index for two
|
binomial proportions using SAS
+--------------------------------------------------------------------------------------------------------------------+
| PARAMETER: N1
: Sample Size 1
|
N2
: Sample Size 2
|
X1
: Outcome 1
|
X2
: Outcome 2
|
Alpha1 : Hyperparameter_α1
|
Alpha2 : Hyperparameter_α2
|
Beta1 : Hyperparameter_β1
|
Beta2 : Hyperparameter_β2
|
Delta : Non-Inferiority Margin
|
Loop : The number of simulations in the group
+--------------------------------------------------------------------------------------------------------------------+
*/
%macro Noninf( n1 = , n2 = , x1 = , x2 = , alpha1 = , alpha2 = , beta1 = , beta2 = , delta = ,
loop = ) ;
%* CREATE FUNCTIONS *;
options cmplib = WORK.FUNC mprint ;
proc fcmp outlib = WORK.FUNC.SAMPLE ;
%* CREATE POCHHAMMER FUNCTION *;
function POCHHAMMER( val , ind ) ;
if ind = 0 then RES = 1 ;
else RES = val ;
do j = 1 to ind - 1 ;
RES = RES * ( val + j ) ;
end ;
return( RES ) ;
endsub ;
%* CREATE FDELTA FUNCTION *;
function FDELTA( i , j , delta , a1 , a2 , b1 , b2 ) ;
do m = 0 to i ;
do n = 0 to j ;
_L1 = POCHHAMMER( 1 - b1 , m ) ;
_L2 = POCHHAMMER( 1 - a2 , n ) ;
_K1 = gamma( a1 + m ) / gamma( a1 ) ; _K2 = gamma( b2 + n ) / gamma( b2 ) ;
_D = ( 1 + delta ) ** m ;
_E = ( 1 + delta ) ** n ;
_U1 = a1 + b2 + m + n ;
_U2 = a1 + b2 ;
_K = ( beta( a1 , b2 ) * ( ( 1 + delta ) ** ( a1 + b2 - 1 ) ) ) /
( beta( a1 , b1 ) * beta( a2 , b2 ) ) ;
if _K1 = 0 or _K2 = 0 or _L1 = 0 or _L2 = 0 or _D = 0 or _E = 0 then KAI = 0 ;
else KAI = ( _K1 * _K2 * _L1 * _L2 * _D * _E ) /
( gamma( _U1 ) / gamma( _U2 ) * fact( n ) * fact( m ) ) ;
ADD + KAI ;
PDF = ADD * _K ;
end ;
end ;
return( PDF ) ;
endsub ;
run ;
%let a1 = %eval( &alpha1 + &x1
);
%let a2 = %eval( &alpha2 + &x2
);
%let b1 = %eval( &n1 - &x1 + &beta1 ) ;
%let b2 = %eval( &n2 - &x2 + &beta2 ) ;
data WORK._FDEL ;
do delta = -0.999 to -&delta by 0.001 ;
PDF = FDELTA( &Loop , &Loop , delta , &a1 , &a2 , &b1 , &b2 ) ;
output ;
end ;
run ;
data WORK.OUT2 ;
set WORK._FDEL end = EOF ;
LAGPDF = LAG( PDF ) ; LAGDEL = LAG( delta ) ;
if _N_ ^= 1 then DAI = sum( LAGPDF , PDF ) * ( delta - LAGDEL ) / 2 ;
RDAI + DAI ;
Tau_e = round( 1 - RDAI , 0.001 ) ;
if EOF ;
run ;
data WORK.OUT3 ;
post1 = &a1 / ( &a1 + &b1 ) ; post2 = &a2 / ( &a2 + &b2 ) ;
deno1 = ( post1 * ( 1 - post1 ) ) / ( &a1 + &b1 + 1 ) ;
deno2 = ( post2 * ( 1 - post2 ) ) / ( &a2 + &b2 + 1 ) ;
random = ( (-1) * &delta - ( post1 - post2 ) ) / sqrt( deno1 + deno2 ) ;
Tau_a = 1 - round( cdf( 'NORMAL' , random , 0 , 1 ) , 0.001 ) ; keep Tau_a ;
run ;
data WORK.OUT4 ;
format n1 n2 x1 x2 alpha1 alpha2 beta1 beta2 delta Tau_a Tau_e ;
merge WORK.OUT2
WORK.OUT3 ;
alpha1 = &alpha1 ; alpha2 = &alpha2 ;
beta1 = &beta1 ; beta2 = &beta2 ;
n1 = &n1 ; n2 = &n2 ; x1 = &x1 ; x2 = &x2 ; a1 = &a1 ; a2 = &a2 ; b1 = &b1 ; b2 = &b2 ;
delta = &delta ;
run ;
proc transpose data = WORK.OUT4 out = WORK.Trs ;
var n1 n2 x1 x2 alpha1 alpha2 beta1 beta2 delta Tau_a Tau_e ;
run ;
data WORK.OUTDS ;
set WORK.Trs ;
length ITEM $ 100 ;
if
_NAME_ = "n1"
then ITEM = "Sample Size 1" ;
else if _NAME_ = "n2"
then ITEM = "Sample Size 2" ;
else if _NAME_ = "x1"
then ITEM = "Outcome 1" ;
else if _NAME_ = "x2"
then ITEM = "Outcome 2" ;
else if _NAME_ = "alpha1" then ITEM = "Hyperparameter_α1" ;
else if _NAME_ = "alpha2" then ITEM = "Hyperparameter_α2" ;
else if _NAME_ = "beta1“ then ITEM = "Hyperparameter_β1" ;
else if _NAME_ = "beta2" then ITEM = "Hyperparameter_β2" ;
else if _NAME_ = "delta" then ITEM = "Non-Inferiority Margin" ;
else if _NAME_ = "Tau_a" then ITEM = "Approximate Probability" ;
else if _NAME_ = "Tau_e" then ITEM = "Exact Probability" ;
if
_NAME_ in( "Tau_a" "Tau_e" ) then VALUE = compress( put( COL1 , 8.3 ) ) ;
else if _NAME_ = "delta"
then VALUE = compress( put( COL1 , 8.2 ) ) ;
else
VALUE = compress( put( COL1 , 8.0 ) ) ;
proc print noobs label ; var ITEM VALUE ; label ITEM = "Item" VALUE = "Value" ; run ;
%mend ;
%Noninf( n1 = 10 , n2 = 10 , x1 = 1 , x2 = 1 , alpha1 = 2 , alpha2 = 3 ,
beta1 = 4 , beta2 = 5 , delta = 0.1 , loop = 100 )
Conclusion
□ 本研究のまとめ
・母比率の差の非劣性を表す指標 τ を紹介した
・事後比率の差の正確なPDFを紹介し、正確確率と近似確率の算出方法を紹介した
・近似確率の算出方法は容易でありプログラムも容易であった。
・正確確率の算出ではFCMPプロシジャを用いることで比較的容易に作成できた。
□ τ を利用する利点
１．直感的に理解がしやすい指標である
２．τ の計算を、SASでプログラミング化することも容易である
３．過去の試験結果の情報を活用して、母比率を比較することが可能
（但し，活用する試験結果は慎重に決めるべきである）
□ 結論
母比率の差の非劣性を示すために、 τ は有用な指標であり、FCMPプロシジャを用いること
で容易に計算ができることが分かった。
29
Reference

Kawasaki Y. and Miyaoka E. (2013). A Bayesian non-inferiority test for two
independent binomial proportions. Pharmaceutical Statisitcs. In press.

Kawasaki Y. and E. Miyaoka E.(2010). A Bayesian inference of P(π1 > π2) for two
proportions. Journal of Biopharmaceutical Statistics. to appear.

TJN-318 Solution Study Group. (1992). Double-blind study on TJN-318 cream and
Bifonazole cream to the patient with the cutaneous mycosis. The Nishinihon
Journal of Dermatology. 54:977–992. (in Japanese).

Pharmaceuticals and Medical Devices Agency. (2009). Follitropin alpha of
assessment report. PMDA: Tokyo. (in Japanese)
30

Download Report