Comments
Description
Transcript
シンポジウム「統計科学の新展開]」報告書
ࢩ࣏ࣥࢪ࣒࢘ࠕ⤫ィ⛉Ꮫࡢ᪂ᒎ㛤@ࠖሗ࿌᭩
ᮏࢩ࣏ࣥࢪ࣒࢘ࡣ⛉Ꮫ◊✲㈝࣭ᇶ┙◊✲(A)ࠕ㠀ᑐ⛠࣭㠀⥺ᙧ⤫ィ⌮ㄽ⤒῭࣭⏕య⛉Ꮫ
ࡢᛂ⏝ࠖ㸦௦⾲㸸㇂ཱྀṇಙ㸦᪩✄⏣Ꮫᇶᖿ⌮ᕤᏛ㒊㸧㸧ࡢຓᡂࡼࡾ㛤ദࡉࢀࡓࠋ
1.㛤ദ᪥㸸2013 ᖺ 11 ᭶ 27 ᪥㸦Ỉ㸧㹼11 ᭶ 29 ᪥㸦㔠㸧
2.ሙ㸸 㔠ἑᏛࢧࢸࣛࢺ࣭ࣉࣛࢨ 3 㝵㞟ᐊ
ୗグࡢ┠ḟࡣㅮ₇㛤ጞ᪥㸪ㅮ₇⪅㸪ㅮ₇⪅ᡤᒓ㸪㡫㸪₇㢟ࡢࡳᥖ㍕ࡋࠊඹྠ◊✲⪅ࡢ
ሗࡣ┬ࡃࠋ
┠ḟ
᭶ ᪥
ᒾୗ Ⓩᚿஓ
ᮾி⌮⛉Ꮫ⌮ᕤᏛ㒊ᩍ㣴
SS
$1(&(66$5<7(67)25(//,37,&$/6<00(75<%$6('2181,)250',675,%87,2129(567,()(/
0$1,)2/'
ᐑ⏣ᗤ୍
㧗ᓮ⤒῭Ꮫ⤒῭Ꮫ㒊
SS
$V\PSWRWLFSURSHUWLHVRI%D\HVLDQW\SHHVWLPDWRUVZLWKRXWDVVXPLQJWKDWWKHKHVVLDQ
PDWULFHVRIORJOLNHOLKRRGIXQFWLRQVFRQYHUJH
ᰠἑె⚽
ᾏ㐨Ꮫ⤒῭Ꮫ㒊
SS
7KLUGRUGHU DYHUDJH ORFDO SRZHUV RI %DUWOHWWW\SH DGMXVWHG WHVWV2UGLQDU\ ^?LW
YHUVXV`DGMXVWHGSURILOHOLNHOLKRRG
⻄ᕝ㞞அ
ᦤ༡Ꮫ⤒῭Ꮫ㒊
SS
7ZR6DPSOH(VWLPDWLRQRI9DU\LQJ&RHIILFLHQW0RGHOVYLD1HDUHVW1HLJKERU0DWFKLQJ
ᒣᮧ⬟㑻
᫂ᏛᏛ㝔ࢢ࣮ࣟࣂ࣭ࣝࣅࢪࢿࢫ◊✲⛉
(PSLULFDO&UHGLW5LVN$QDO\VLVRQ(XUR*RYHUQPHQW%RQGV
I
SS
ᒣ᪉ு
ୖᬛᏛᏛ㝔⌮ᕤᏛ◊✲⛉ᩘᏛ㡿ᇦ
SS
2SWLPDOLQYHVWPHQWLQFRUUHODWHGVWRFNVDQGDQLQGH[ERQGIRUDGHILQHGFRQWULEXWLRQ
SHQVLRQ
㛗⏫▱᠇
ୖᬛᏛᏛ㝔⌮ᕤᏛ◊✲⛉ᩘᏛ㡿ᇦ
SS
ࢧ࣮࢟ࢵࢺࣈ࣮࣮ࣞ࢝ไᗘࢆ⪃៖ࡋࡓඛ≀ၟရࡢࣜࢫࢡホ౯
1JDL+DQJ
&KDQ
&KLQHVH8QLYRI+RQJ.RQJ
S
6WDWLVWLFDO$UELWUDJHRQ9RODWLOLW\
5RJHU.RHQNHU
8QLY,OOLQRLV S
8QREVHUYHG+HWHURJHQHLW\LQ/RQJLWXGLQDO'DWD$Q(PSLULFDO%D\HV3HUVSHFWLYH
᭶ ᪥
ᙇ ඖ᐀
┠ⓑᏛ♫ሗᏛ⛉
SS
☜⋡ඃ㉺ᛶࡢホ౯ᇶ‽࡛ᖹᆒศᩓࡢ᪉㡰ᗎไ⣙ࡀ࠶ࡿሙྜࡢ㸰ࡘࡢṇつẕ㞟ᅋࡢᖹ
ᆒࡢ᥎ᐃࡑࡢᛂ⏝
ΎỈ⚈ኴ
ᓥ᰿ᏛᏛ㝔⥲ྜ⌮ᕤᏛ◊✲⛉ SS
ᡂ㋣Ꮫ⤒῭Ꮫ㒊
/065PHWKRG ࡑࡢᛂ⏝
᪂ᮧ⚽୍
SS
SS
SS
SS
ศᩓඹศᩓ⾜ิᇶ࡙ࡃุู㛵ᩘࡢ⤊↉
ᒣᮏࡅ࠸Ꮚ
ภ㤋ᕤᴗ㧗➼ᑓ㛛Ꮫᰯ
࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࢆ⏝࠸ࡓ㠀⥺ᙧุูᡭἲࡢᥦ
ゝ ᪩✄⏣Ꮫᇶᖿ⌮ᕤᏛ◊✲⛉
$V\PSWRWLFVIRU0(VWLPDWRUVLQ7LPH6HULHV
⏣୰ே
Ꮫ⩦㝔Ꮫ⤒῭Ꮫ㒊
6WDWLVWLFDOLQIHUHQFHDVVRFLDWHGZLWKWKHIUDFWLRQDO%URZQLDQPRWLRQDQGUHODWHG
SURFHVVHV
II
బ⸨୍
ᗈᓥᏛཎ⇿ᨺᑕ⥺་⛉Ꮫ◊✲ᡤ SS
㛫ࡸኚᩘ✵㛫ୖ࡛ኚࡍࡿᅇᖐಀᩘࡘ࠸࡚
⏣ဴ
┴❧ᗈᓥᏛ⤒ႠሗᏛ㒊
SS
SS
᭷㝈༊㛫࠾ࡅࡿኚಀᩘࡢྠಙ㢗༊㛫ࡢᵓ⠏ࡘ࠸࡚
Ἠᚿὠᜨ
ศᏛᕤᏛ㒊 ⤒ࢸ࢟ࢫࢺࢹ࣮ࢱᑐࡍࡿከḟඖᑻᗘἲࡢᛂ⏝
⳹ᒣᐉ⬍
ᑦ⨾ᏛᅬᏛⱁ⾡ሗᏛ㒊ሗ⾲⌧Ꮫ⛉
SS
ேཱྀືែ⤫ィᇶ࡙ࡃே㛫ࡢᑑ㝈⏺ࡢ᥎ᐃ
᭶ ᪥
㧗ᶫ ᑗᐅ
⊂❧⾜ᨻἲே⤫ィࢭࣥࢱ࣮⤫ィᢏ⾡◊✲ㄢ SS
つᶍ⤒῭⣔ࢹ࣮ࢱ࠾ࡅࡿᵝࠎ࡞ከ㔜௦ධἲࣝࢦࣜࢬ࣒ࡢ᳨ド
ఀ⸨ఙ
᫂ᾏᏛ⤒῭Ꮫ㒊
SS
SS
SS
SS
ᅜໃㄪᰝ࣑ࢡࣟࢹ࣮ࢱࢆ⏝࠸ࡓ༏ྡᢏἲࡢ᭷ຠᛶࡢ᳨ド
ඖ
㮵ඣᓥᏛ
ࣆࢵࢺ࣐ࣥ☜⋡ศࡢศᩘࡢᴟ㝈ศᕸࡢ୍⪃ᐹ
㇂ᨻ
᠕⩏ሿᏛ⌮ᕤᏛ㒊
☜⋡ศࡢᶆᮏண 㔞㸸⏕ែᏛࡢᛂ⏝
ሷ℈ᩗஅ
ᮾி⌮⛉ᏛᕤᏛ㒊
6HPLSDUDPHWULF(IILFLHQF\IRUWKH4XDQWLOH5HJUHVVLRQEDVHG/HVWLPDWLRQ
㇂ྜᘯ⾜
᪩✄⏣Ꮫᅜ㝿ᩍ㣴Ꮫ㒊 SS
2QDWDQJHQWVSDFHIRUWKHFRHIILFLHQWIXQFWLRQVRI4XDQWLOH5HJUHVVLRQ㸦ศᩘ
ᅇᖐࣔࢹࣝࡢࣛࣥࢲ࣒ಀᩘ⾲⌧㛵ࡍࡿ᥋✵㛫ࡘ࠸࡚㸧
⏕டΎ㈗
ᮾி⌮⛉Ꮫ⌮ᕤᏛ㒊
㡰ᗎ࢝ࢸࢦࣜṇ᪉ศ⾲࠾ࡅࡿ㸰ኚ㔞 W ศᕸᆺᑐ⛠ࣔࢹࣝࡘ࠸࡚
III
SS
ᓥ⏣ᩥ㤶
ᮾி⌮⛉ᏛᏛ㝔⌮ᕤᏛ◊✲⛉ SS
㡰ᗎ࢝ࢸࢦࣜṇ᪉ศ⾲࠾ࡅࡿᑐ⛠ᛶࡢࣔࢹࣝศゎ࠾ࡼࡧేྜࡋࡓ⾲ᇶ࡙ࡃᑐ⛠ᛶ
㛵ࡍࡿᑻᗘ
Ᏻ⸨᐀ྖ
ᮾி⌮⛉ᏛᏛ㝔⌮ᕤᏛ◊✲⛉ SS
ṇ᪉ศ⾲࠾ࡅࡿ⣼✚☜⋡ࢆ⏝࠸ࡓ㠀ᑐ⛠ᛶࡢࣔࢹࣝࡢศゎ
୕ᯞ♸㍜
ᮾி⌮⛉ᏛᏛ㝔⌮ᕤᏛ◊✲⛉ ṇ᪉ศ⾲࠾ࡅࡿᣑᙇࣃࣜࣥࢻ࣑ࣟࢵࢡᑐ⛠ࣔࢹࣝᑐ⛠ᛶࡢศゎ
IV
SS
A necessary test for elliptical symmetry based on uniform distribution over Stiefel manifold
·
Karlsruher Institut für Technologie
Bernhard Klar
1.
( )
,
,
.
,
,
,
,
,
X 1, . . . , X N
p-
( )
.
p×N
. X = [X 1 , . . . , X N ]
,
X̄ = N −1 Xj N ,
(1)
S = n−1 XQN X " ,
.
n = N − 1 ≥ p,
(2)
,
QN = IN − N −1 j N j "N .
j N = (1, . . . , 1)" ∈ RN ,
,
(3)
p×N
,
W = [W 1 , . . . , W N ] = S −1/2 [X 1 , . . . , X N ] QN = S −1/2 XQN ,
,
.
W i = S −1/2 (X i − X̄),
S −1/2
(4)
,
i = 1, . . . , N,
(5)
S
2.
X 1, . . . , X N
matrix)
.
X ∼ ECp (0, Λ)
, X = [X 1 , . . . , X N ]
LSp×N (φ)
p×N
(observation
Fp×N
Fp×N = {X(p × N ) ∼ LSp×N (φX ); the distribution of Xa depends on a ∈ RN only through a" a}
1
(6)
Iwashita and Klar (2013)
Theorem 1.
.
X 1, . . . , X N
ECp (0, Λ)
X = [X 1 , . . . , X N ]
Y = S −1/2 X
p×N
. S
(2)
.
Y ∼ SSp×N (φY ).
, (3)
QN
N ×n ( )
,S
, KK " = QN , K " K = In
rank(QN ) = n
K
.
Theorem 1
Theorem 2. Let X 1 , X 2 , . . . , X N
.
(7)
,
ECp (µ, Λ)
.
, X = [X 1 , X 2 , . . . , X N ]
(2)
.
, n×p
U = K " QN X " (nS)−1/2 = K " X " (nS)−1/2
O(n, p)
Stiefel
.
, QN
(8)
N ×N
(3)
.
3.
(k)
{X i }N
i=1 (k = 1, 2, . . . , m)
X(k)
p!
"
(k)
(k)
= X1 , . . . , XN ,
X
"
S(k) = n−1 X(k) QN X(k)
,
"
"
Uk = K(k)
X(k)
(nS(k) )−1/2 ,
(k)
, Theorem 2
, Xj
,
n = N − 1 ≥ p,
∼ ECp (µ, Λ) (j = 1, . . . , N ; k = 1, . . . , m)
O(n, p)
Jupp (2001)
,
,
.
Uk
Stiefel
Stiefel
,
,
,
(9)
,
.
.
Iwashita, T. and Klar, B.(2013). The joint distribution of Studentized residuals under elliptical distributions.
submitted
Jupp, P. E.(2001). Modification of the Rayleigh and Bingham tests for uniformity of directions. J. Multivariate Anal., 77, 1–20.
2
Asymptotic properties of Bayesian type estimators without
assuming that the Hessian matrices of log-likelihood functions
converge
(Van der Vaart 1998,
Ibragimov and Hasminskii 1981, Ferguson 1996, Basawa and Rao 1980, Yoshida 2007)
.
,
Bernstein-von Mises theorem
,
,
.
White (1980, pp.727-728)
!
.
n
Θ ≡ (Θ̃1 , ..., Θ̃d )!
Θ̃
Ỹ
:
Qn (Ỹ, θ )
Qn (θθ )
b(θθ )
(1)
pn (θθ |Ỹ) = !
exp{−nQn (Ỹ, θ )}b(θθ )
.
θ
θ
exp{−nQ
(
Ỹ,
θ
)}b(θ
)dθ
n
Θ
Ln (θθ )
θ
, π(θθ )
Qn (θθ )
−(1/n) log Ln (θθ )π(θθ )
Θ
Θ̃
.
Ỹ
(1)
!
"!
#!
θ1 exp{−nQn (θθ )}b(θθ )dθθ
θk exp{−nQn (θθ )}b(θθ )dθθ
Θ
Θ
!
≡
, ..., !
.
θ
θ θ
θ
θ θ
Θ exp{−nQn (θ )}b(θ )dθ
Θ exp{−nQn (θ )}b(θ )dθ
θB
n
Hanousek (1990)
,
, Qn (θθ ) = −(1/n) log Ln (θθ ), Qn (θθ ) =
.
.
.
(1)
i.i.d.
(2)
Qn (θθ )
,
.
, Qn (θθ )
.
,
,
0
.
.
α-mixing
Gourieroux and Monfort (1981)
.
AR(1)
, α-mixing
.
,
AR(1)
α-mixing
.
Third-order average local powers of Bartlett-type adjusted tests:
Ordinary versus adjusted profile likelihood
柿沢 佳秀 (北大経済)
1. 設定と記号 dX -次元の密度関数モデル P = {f (x, θ) : θ ∈ Θ ⊂ Rp } で未知の
母数ベクトルを θ = (θ !(1) , θ !(2) )! , θ (1) = (θ1 , . . . , θp1 )! , θ (2) = (θp1 +1 , . . . , θp )! のように
分割し, 帰無仮説 θ (1) = θ (1)0 を検定する. "ji1 ...jiRi (x, θ) = ∂θ∂j · · · ∂θj∂ log f (x, θ),
i1
iRi
i = 1, . . . , v の f (x, θ) に関する v 次キュムラントを νj11 ...j1R1 ,...,jv1...jvRv (θ) と書き,
!
(N )
!
(N )
N
−1/2
Zj (θ) = N −1/2 N
i=1 "j (Xi , θ), Zj1 j2 (θ) = N
i=1 {"j1 j2 (Xi , θ) − νj1 j2 (θ)} などを
定義する. このとき, 情報量行列 [νj,k (θ)] = −[νjk (θ)] とスコアベクトルを分割し,
[νj,k (θ)] =
"
ν (1,1) (θ) ν (1,2) (θ)
ν (2,1) (θ) ν (2,2) (θ)
#
,
(N )
[Zj (θ)]
(N )
Z (θ)
= (1)
(N )
Z(2) (θ)
と書く. ここで, パラメータ直交性, すなわち
, 『任意の θ ∈ Θ
に対し ν (1,2) (θ) = Op1 ,p2 』
Ip 1
を定義する. 以下, 添え字
を仮定しない. p×p1 行列 [Gj,a (θ)] =
−1
−ν (2,2) (θ)ν (2,1) (θ)
については, j1 , . . . , jR ∈ {1, . . . , p} (j, k も同様), a1 , . . . , aR ∈ {1, . . . , p1 } (a, b も同様),
!
!
r1 , . . . , rR ∈ {p1 + 1, . . . , p} (r, s も同様) と約束し, それらの和 pj1 ,...,jR =1 , pa11 ,...,aR =1 ,
!p
r1 ,...,rR =p1 +1 に対しては Einstein 表記法を採用する.
2. 通常の尤度解析 対数尤度 L(N ) (θ) =
MLE
(N )
θ(
及び Restricted
LR(N ) =
R(N ) =
!N
i=1
log f (Xi, θ) を最大化し, (Unrestricted)
(N )
MLE θ)
での評価は, ハット/チルド記号を付ける;
( (N ) − θ
) (N ) )! N 1/2 (θ
2(L((N ) − L)(N ) ) , grad(N ) = (Z
)
(1)0
(1)
(1)
(N
)
(N
)
d
(N
)
(N
)
! ) −1
(N )
!(
(
(
)
)
(Z(1) ) ν (11·2) Z(1) , W
= N(θ (1) − θ (1)0 ) ν (11·2) (θ (1) − θ (1)0 ) −→
(N )
!)
( (N )
( (N )
) (N ) ! ( −1 ) (N )
MR(N ) = (Z(1) ) ν (11·2) Z(1) , MW
χ2p1 .
= N(θ (1) − θ (1)0 ) ν (11·2) (θ (1) − θ (1)0 )
3. 修正プロファイル尤度解析 各 θ (1) に対しプロファイル対数尤度は
LP(N ) (θ (1) ) = L(N ) (θ̌
で定義される. ここに, θ̌
(N )
(θ (1) ))
θ (1)
, θ̌ (N ) (θ (1) ) = arg maxθ L(N ) (θ)
(N )
(2)
(2)
θ̌ (2) (θ (1) )
(N )
θ̌ (θ (1)0 ) に注意すれば, 上記の LR(N ) は
(θ (1) ) =
( (N ) = θ̌ (N ) (θ
( (N ) ) と θ
) (N ) =
であり, θ
(1)
(N )
(N )
d
LR(N ) = 2{LP(N ) (θ( (1) ) − LP(N ) (θ (1)0 )} −→ χ2p1
と書け, 局外母数を含む場合の尤度比検定は “プロファイル尤度解析”と解釈されうる.
(N )
しかし, LP(N ) (θ (1) ) は “正真正銘な対数尤度”でなく, 『Eθ ,θ† [ ∂θ∂a LP(N ) (θ (1)0 )] = 0』
(1)0
(2)
とならないため, Stern (1997; JRSS) は修正プロファイル対数尤度を定義している;
(N )
LAP(N ) (θ (1) ) = LP(N ) (θ (1) ) − N −1/2 M̌! Ǧ ν̌ −1
(11·2) Ž(1)
[ここに, Q̌ = Q(θ̌
(N )
(θ (1) ))].
r,r !
なお, Mj (θ) = 12 {νrr! ,j (θ) + νr,r! ,j (θ)}ν(2,2) (θ) は局外母数を含む高次漸近論の文脈で
お馴染みのものであり, この量が消えないとき LR(N ) は一般的に２次漸近不偏でない
(Hayakawa, 1975; Bio). 本報告では, 2 節の通常の統計量とその修正プロファイル版
(N )
LRAP(N ) = 2{LAP(N ) (θ (1) ) − LAP(N ) (θ (1)0 )} ,
(N )
gradAP(N ) = [ ∂θ∂a LAP(N ) (θ (1)0 )][θ (1) − θ (1)0 ]a ,
a1 ,a2
RAP(N ) = N −1 [ ∂θ∂a LAP(N ) (θ (1)0 )]ν(11·2)
(θ̌
(N )
1
W
AP(N )
(N )
= N(θ (1)
!
− θ (1)0 ) ν (11·2) (θ̌
(N )
2
(N )
(N )
(θ (1) ))(θ (1)
a1 ,a2
(θ̌
MRAP(N ) = N −1 [ ∂θ∂a LAP(N ) (θ (1)0 )]ν(11·2)
(N )
1
MW
AP(N )
(N )
= N(θ (1)
− θ (1)0 )! ν (11·2) (θ̌
(N )
(θ (1)0 ))[ ∂θ∂a LAP(N ) (θ (1)0 )] ,
− θ (1)0 ) ,
(N )
(θ (1) ))[ ∂θ∂a LAP(N ) (θ (1)0 )] ,
2
(N )
(θ (1)0 ))(θ (1) − θ (1)0 )
(N )
の局所検出力を漸近展開で比較した. ここに, θ (1) = arg maxθ(1) LAP(N ) (θ (1) ).
4. (平均) 局所検出力 2,3 節の統計量は, それぞれ
)(N ) +
)a(N ) ν)a,b Z
T (N ) ≈ Z
(11·2) b
N
3
2
0
/
2 . ) G G G / −1 ) (N )
(N )
GG
) (N )
)
)
Z
C
[)
ν
]
+
C
[)
ν −1
Z
(11·2) (1) bi
(11·2) Z(1) ]bi
b1 b2 b3
b1 b2 ,k1 k2 k1 k2
1/2
i=1
i=1
4
3
/
2 1 ) G G G G / −1 ) (N )
) (N )
)(N ) + D
)(N ) ) [)
)G G G
)G G G
+
Z
Z
Db1 b2 b3 b4
[)
ν (11·2) Z(1) ]bi + (D
ν −1
(11·2) Z(1) ]bi
b1 b2 b3 ,k1 k2 k1 k2
b1 b2 b3 ,k1 k2 k3 k1 k2 k3
N
i=1
i=1
(N )
T
AP(N )
≈T
+
(N )
)
)
)G G
+D
b1 b2 ,k1 k2 ,k3 k4 Zk1 k2 Zk3 k4
(N )
+
2
2
/
i=1
3G [)
) (N ) ]b
M
ν −1 Z
N 1/2 b1 (11·2) (1) 1
2
) (N )
[)
ν −1
(11·2) Z(1) ]bi ,
2
0
/
2 . 1 3G b1 ,b2 3G
)G G
) (N )
)G
)(N ) ν −1 Z
) (N )
Mb1 ν)(11·2) Mb2 + M D
[)
ν −1
(11·2) Z(1) ]bi + M Db1 ,k1 k2 Zk1 k2 [)
(11·2) (1) ]b1
b1 b2
N 2
i=1
と確率展開され, 以下の２形式の調整 (Kakizawa, 2012ab,2013; JMVA,SPL) を考えた;
#GCF(N ) = #(N ) +
#GB(N ) = #(N ) +
2
N
4
R=2,4,6
)
)
+∆
b1 b2 b3 ,k1 k2 Zk1 k2
(N )
†
(1)0 ,θ(2)
i=1
(N )
)
[ν) −1
(11·2) Z(1) ]bi ,
4
1
/
2 )C
−1
) (N ) ] + 2 (Γ
)#
)
) (N )
)
Γ
[
ν
Z
+
∆
)
[ν) −1
b
b1 b2 b3 b4
b1 b2 b3 b4
(11·2) Z(1) ]bi
N 1/2 b1 b2 b3 i=1 (11·2) (1) i N
i=1
(N )
ただし, Pθ
R
/
)
Γ
b1 ···bR
3
/
3
/
(N )
)
)#
[ν) −1
(11·2) Z(1) ]bi + Γb1 b2
i=1
(N )
[#GCF(N ) ≤ x], Pθ
†
(1)0 ,θ(2)
2
/
(N )
2
)
[ν) −1
(11·2) Z(1) ]bi .
i=1
[#GB(N ) ≤ x] = Pr[χ2p1 ≤ x] + o(N −1 ) となる
ようにした (こうしないとサイズが揃わないために公平な検出力比較が不可能である).
局所対立仮説 θ †(N ) ≡ (θ !(1)0 , (θ †(2) )! )! + N −1/2 (h!(1) , 0!p2 )! の下で検出力を 3 次まで導出
して, 以下の平均基準から RGCF(N ) または (RAP )GCF(N ) の “ある意味”の 3 次最適性を
示した;
ave{π(h(1) )} =
Sλ
5
Sλ
π(h(1) ) dh(1)
5
,
Sλ dh(1)
Sλ = {h(1) ∈ Rp1 : h!(1) ν (11·2) h(1) = λ} (λ > 0) .
Ặ
ྡ㸸⻄ᕝ㞞அ
ᡤ
ᒓ㸸ᦤ༡Ꮫ⤒῭Ꮫ㒊
ㅮ₇㢟┠ 㸸 “Two-Sample Estimation of Varying Coefficient Models via Nearest
Neighbor Matching” [Artem Prokhorov (University of Sydney)ࡢඹⴭ]
ㅮ₇ෆᐜ㸸
᭱ึ᩿ࡗ࡚࠾ࡃࠊㅮ₇㢟┠ࡣྍኚಀᩘᅇᖐࣔࢹࣝ(varying coefficient
model, “VCM”)ࡢ᥎ᐃ⌮ㄽ㛵ࡍࡿࡶࡢ࡛࠶ࡿࡀࠊ」ᩘࡢࢹ࣮ࢱࢭࢵࢺࡽ᭱
㏆㞄࣐ࢵࢳࣥࢢ(nearest neighbor matching, “NNM”)ࢆᛂ⏝ࡋ࡚᪂ࡓࢹ࣮ࢱࢭࢵ
ࢺࢆᵓ⠏ࡍࡿሙྜࡢ᥎ᐃ⌮ㄽࡣ୍⯡࡞ࡌࡳࡀⷧ࠸ணࡉࢀࡓࠋࡑࡢࡓࡵࠊ
ᮏㅮ₇࡛ࡣࠊ⥺ᙧᅇᖐࣔࢹࣝ㛵ࡍࡿ᥎ᐃ⌮ㄽࢆヲ⣽ㄝ᫂ࡋࠊࡑࡢ⌮ㄽࡢᣑ
ᙇࡋ࡚ࠊVCM 㛵ࡍࡿ᥎ᐃ⌮ㄽࢆᴫㄝࡍࡿ࠸࠺ᡭ㡰ࢆ᥇⏝ࡋࡓࠋ
ᅇᖐࣔࢹࣝ㛵ࡍࡿ᥎ᐃ⌮ㄽࡣࠊࣔࢹࣝྵࡲࢀࡿኚᩘ࡚ࡀ༢୍ࡢࢹ࣮ࢱ
ࢭࢵࢺࡽᚓࡽࢀࡿࡇࢆ๓ᥦࡋ࡚ᵓ⠏ࡉࢀ࡚࠸ࡿࠋࡋࡋࠊປാ⤒῭Ꮫ࣭
බඹ⤒῭Ꮫ࡞ࡢᐇドศᯒ࡛ࡣࠊᅇᖐࣔࢹࣝࡢ᥎ᐃᚲせ࡞ࡿኚᩘࡀ」ᩘࡢ
ࢹ࣮ࢱࢭࢵࢺࢃࡓࡗ࡚Ꮡᅾࡍࡿࢣ࣮ࢫࡀࡋࡤࡋࡤ࠶ࡾࠊࡇࡢሙྜࠊࡇࢀࡽࡢ
ࢹ࣮ࢱࢭࢵࢺࡽ NNM ࢆ⏝࠸࡚ศᯒᚲせ࡞ࢹ࣮ࢱࢭࢵࢺࢆసᡂࡍࡿࡇࡀ
័⾜ࡋ࡚⾜ࢃࢀ࡚࠸ࡿࠋࡇࡢ័⾜ࡣࡉࡽࠊᅇᖐࣔࢹࣝྑ㎶ࡢㄝ᫂ኚᩘయ
ࡀ NNM ࡛࣐ࢵࢳࡉࡏࡓࢹ࣮ࢱ⨨ࡁ࠼ࡽࢀࡿࢣ࣮ࢫࠊㄝ᫂ኚᩘࡢ୍㒊ࡢ
ࡳࡀ⨨ࡁ࠼ࡽࢀࡿࢣ࣮ࢫࡢ㸰㏻ࡾศ㢮ྍ⬟࡛࠶ࡿࠋศᯒࡢᐜ᫆ࡉࡽࠊ
ᅇࡣ๓⪅ࡢࡳࢆ⪃ᐹࡢᑐ㇟ࡋࡓࠋ
ᮏㅮ₇࡛ࡣࠊࡲࡎࠊNNM ࡛ᵓ⠏ࡋࡓࢹ࣮ࢱࢭࢵࢺࢆ⏝࠸࡚⥺ᙧᅇᖐࣔࢹࣝ
ࡢಀᩘࢆ᭱ᑠἲࡼࡾ᥎ᐃࡍࡿሙྜࠊ᭱ᑠ(two-sample ordinary least
squares, “TSOLS”)᥎ᐃ㔞ࡣ୍⮴ᛶࢆᣢࡓ࡞࠸ࡇࢆ♧ࡋࡓࠋ┤ឤⓗࡣࠊ⿕ㄝ
᫂ኚᩘㄝ᫂ኚᩘࡢ࣑ࢫ࣐ࢵࢳ㉳ᅉࡍࡿࣔࢹࣝᐃᘧࡢㄗࡾࡢ⤖ᯝࠊ㝖እ
ኚᩘࣂࢫ(omitted variable bias)㢮ఝࡍࡿࣂࢫࡀⓎ⏕ࡋࠊTSOLS ᥎ᐃ㔞
ࡣ୍⮴࡞ࡿㄝ᫂ࡉࢀࡿࠋ
TSOLS ᥎ᐃ㔞ࡣ୍⮴ᛶࢆᣢࡓ࡞࠸༙㠃ࠊࡑࡢ☜⋡ᴟ㝈ࡣ┿್ࡢ୍ḟኚࡋ
࡚⾲⌧ࡉࢀࡿࠋࡇࡢ୍ḟኚᑐࡍࡿ㏫ീࡽᵓ⠏ࡉࢀࡿ᥎ᐃ㔞ࡣ┿್ᑐ
1
ࡋ୍࡚⮴ᛶࢆᣢࡘࠋࡇࡢⓎࡣࠊGouriéroux, Monfort and Renault (1993)࠾ࡼࡧ
Smith (1993)ࡢ indirect inference ᥎ᐃἲࡢ୍࡛࠶ࡿࠋࡑࡢࡓࡵࠊࡇࡇ࡛ᥦࡉ
ࢀࡿࣂࢫಟṇ᥎ᐃ㔞ࢆ two-sample indirect inference (“TSII”)᥎ᐃ㔞ࡪࡇ
ࡋࡓࠋTSII ᥎ᐃ㔞ࡢ୍⮴ᛶࠊ㏆ṇつᛶ࠾ࡼࡧ㏆ศᩓࡶྠ♧ࡉࢀࡓࠋ
ࡉࡽࠊࣔࣥࢸ࢝ࣝࣟᐇ㦂ࡼࡾࠊTSOLS ᥎ᐃ㔞ࡢ୍⮴ᛶ࠾ࡼࡧ TSII ᥎ᐃ
㔞ࡢ୍⮴ᛶࡀᩘ್ⓗ☜ㄆࡉࢀࡓࠋTSII ᥎ᐃ㔞ࡣ୍⮴ᛶࢆᣢࡘ㠃ࠊ࡚ࡢኚ
ᩘࡀ༢୍ࡢࢹ࣮ࢱࢭࢵࢺࡽᚓࡽࢀࡿ௬ᐃࡋࡓሙྜࡢᐇ⾜ྍ⬟࡞ OLS ᥎ᐃ
㔞ẚ࡚᭷ຠᛶࡣຎᚋࡍࡿࠋࡇࢀࡣ」ᩘࡢࢹ࣮ࢱࢭࢵࢺࢆ⏝ࡍࡿࡇక
࠺ሗࡢࣟࢫᑐࡍࡿࢥࢫࢺゎ㔘ࡉࢀࡿࠋ
ࡇࡢ⤖ᯝࢆⓎᒎࡉࡏࠊNNM ࡛ᵓ⠏ࡋࡓࢹ࣮ࢱࢭࢵࢺࢆ⏝࠸࡚ VCM ࢆ࣮࢝ࢿ
ࣝᖹࡼࡾ᥎ᐃࡍࡿሙྜࡢ⌮ㄽࡶᚓࡽࢀࡓࠋලయⓗࡣࠊ⥺ᙧᅇᖐࣔࢹࣝ
ࡢሙྜྠᵝࠊྍኚಀᩘࡢᒁᡤ⥺ᙧ(local linear, “LL”)᥎ᐃ㔞ࡣྛⅬ୍࡛⮴ᛶࢆ
ᣢࡓ࡞࠸ࡇࡀド᫂ࡉࢀࡓࠋࡋࡋࠊࡑࡢ☜⋡ᴟ㝈ࡀࡸࡣࡾ┿್ࡢ୍ḟኚ
ࡋ࡚⾲⌧ࡉࢀࡿࡇᇶ࡙ࡁࠊࣂࢫࢆಟṇࡋࡓࣀࣥࣃ࣓ࣛࢺࣜࢵࢡ᥎ᐃ㔞
ࢆᗄࡘࡢࣂ࣮ࢪ࡛ࣙࣥᥦࡋࠊࡇࢀࡽࢆ⥲⛠ࡋ࡚ bias-corrected LL (“BCLL”)
᥎ᐃ㔞ࡪࡇࡋࡓࠋࡇࢀࡽ BCLL ᥎ᐃ㔞㛵ࡋࠊ୍⮴ᛶࠊ㏆ṇつᛶࠊ
ࣂࢫ㡯ࡢ㏆ఝᘧࠊ࠾ࡼࡧ㏆ศᩓࡶ᫂♧ࡉࢀࡓࠋ
᭱ᚋࠊᚋࡢ◊✲ࡢ᪉ྥࡀࠊ(i) BCLL ᥎ᐃ㔞ࡢᐇ⾜ἲ㸦ࣂࣥࢻᖜ㑅ᢥἲ㸧
࠾ࡼࡧࣔࣥࢸ࢝ࣝࣟᐇ㦂ࠊ(ii) ᅇᖐࣔࢹࣝྑ㎶ࡢㄝ᫂ኚᩘࡢ୍㒊ࡢࡳࡀ NNM ࡛
࣐ࢵࢳࡉࡏࡓࢹ࣮ࢱ⨨ࡁ࠼ࡽࢀࡿࢣ࣮ࢫᑐࡍࡿ᥎ᐃ⌮ㄽࡢᣑᙇࠊࡉࡽ
ࡣࠊ(iii) ᐇ㝿ࡢࢹ࣮ࢱࢆ⏝࠸ࡓᐇドศᯒ⛣⾜ࡍࡿணᐃ࡛࠶ࡿࡇゐࢀࠊᮏ
ㅮ₇ࢆ⥾ࡵࡃࡃࡗࡓࠋ≉ࠊ(ii)ࡣࠊࢧ࣮࣭࣋ࢧࣥࣉࣜࣥࢢ࠾࠸࡚ኚᩘࡢㄪ
ᰝ₃ࢀࡀⓎ⏕ࡋࡓሙྜࠊࡢࡼ࠺㏣ຍㄪᰝࢆ⾜࠼ࡤࡼ࠸ࢆ♧၀ࡍࡿෆᐜ
ࡶゎ㔘࡛ࡁࡿࡓࡵࠊࡇࡢ᪉㠃ࡢ◊✲ᣑᙇὀຊࡋࡓ࠸ࠋ
ཧ⪃ᩥ⊩㸸
Gouriéroux, C., A. Monfort, and E. Renault (1993): “Indirect Inference,” Journal of
Applied Econometrics, 8, S85 - S118.
Smith, A. A., Jr. (1993): “Estimating Nonlinear Time-Series Models Using Simulated
Vector Autoregressions,” Journal of Applied Econometrics, 8, S63 - S84.
2
E mpirical C redit Risk Analysis on E uro Government Bonds
± Term Structures of Default Probabilities±
Takeaki Kariya (GSGB, Meiji University) [email protected]
Yoshiro Yamamura(GSGB, Meiji University) [email protected]
Yoko Tanokura (AMS, Meiji University) [email protected]
Zhu Wang（ZW System）
The Financial Crisis in 2008-2009 and the European Crisis thereafter made many
European states (countries) with currency integration or equivalently member states of the
European Economic and Monetary Union (EEMU) confronting severe budgetary and
unemployment problems, as reflected in Greek economy, where as of May 2013, there are
27 countries in the EU (European Union) among which 17 states form the EEMU. The
European problem is affecting global economies through trade relations and financial
markets, and naturally the world concerns about their future movements, because a collapse
of the EEMU would make a significantly serious impact on the world economy.
In this paper, via interest rate (IR) differential, GB price differential, default probability
(DP) and CDS, we will make a comprehensive credit risk analysis on the GB price data of
the Five States; Germany, France, Italy, Spain and Greece over 2007.4-2012.3. In order to
derive term structures of interest rates (TSIRs), the forward-looking GB-pricing model
proposed in Kariya (1993) and applied in Kariya, Wang, Wang, Doi and Yamamura (2012)
(shortened KWWDY (2012)) is applied to each monthly cross-sectional set of GB prices.
And to derive term structures of default probabilities (TSDPs) relatively to DGB, the
forward-looking CB (Corporate Bond) pricing model proposed in Kariya (2013) is also
applied to each monthly cross-sectional set of FGB, IGB, SGB, and GrGB, which are
regarded as CBs in the model, while DCBs are viewed as non-defaultable reference GBs. In
derivations of TSIRs, P-differentials and TSDPs, we use each monthly cross-sectional set
of GB price data at the last business day of each month where the period of analysis is
2007.4-2012.3. Our arguments are sometimes made associated with the conditions in
Maastricht Treaty, which we will discuss in Section 2 and show a legitimacy of regarding
DGB as a reference GB.
First after making some discussions and observations on the business cycles of the
Five States and the Maastricht convergence condition of IRs, we derived TSIRs via the
bond-pricing model in KWWDY (2012) and compared them. In association with the
budgetary conditions and business cycles of the Five States, we found that (a) in the
sub-period 2007.4-2008.6 the TSIRs of the Five States moved together at the
approximately same levels, meaning that no credit differentiation was found in
IR-differentials, (b) in the sub-period 2008.7-2010.3 all the IRs tended to decrease
gradually as trend movements but their IR-differentials exhibited some credit
differentiation, (c) in the sub-period 2010.4-2011.12 the IR-differentials tended to diverge
and (d) after that they tended to get stabilized gradually. Though F-TSIRs moved together
with D-TSIR before the Euro Crisis, then they deviate. This deviation would have made the
role of the German Government more important for the stability of the EEMU.
Secondly we proposed what we call CRPS measure, which measures directly credit risk
via P-differential between a given GB and DGB-equivalent in terms of euros. The CRPS
measure is model-free once the mean discount function of Germany is estimated. In our
empirical analysis the 10-year CRPSs of Italy and Spain in the mid of the Financial Crisis
2009.3 were respectively about -10 euros and -7 euros, but in the mid of the Euro Crisis
2011.12 they jumped down to about -35 euros and -30 euros. Since the CRPS measure is
additive, it can be used to measure credit risk volume (in euros) of a bond portfolio.
Thirdly the TSDPs of the Four States were derived via CB-pricing model in Kaiya (2013)
and compared. The model enabled us to transform the term structures of CRPSs into
TSDPs, and it turned out that the 10 year CRPSs of Italy and Spain in 2011.12 respectively
corresponded to about 35% and 30% where the recovery rate was assumed to be zero. In
addition we substantiated the observations on credit differentiation obtained through TSIR
analysis in terms of DPs. Also it was observed that the Financial Crisis did not affect the
DPs of France, did increase the DPs of Italy, Spain and Greece to the levels of 6%, 10%
and 20% respectively, but after that their DPs decreased. The time series movements of the
TSDPs of the Four States were associated with business cycles, Financial Crisis and Euro
Crisis.
Fourthly the time series relationships between IR-differentials and DPs of each maturity
were shown to be strongly linear by our regression analysis, which enables us to convert
IR-differentials into DPs for each maturity and vice versa. Consequently the Maastricht
convergence condition can be stated in terms of DPs. It was observed that Italy, Spain and
Greece did not meet the required condition in the period of the Euro Crisis though they are
the members of the EEMU, creating the instability of the EEMU system that was an
economic concern in the global world. But an explicit solution that President Van Rompuy
(2012) planned for a genuine integration of the EEMU will stabilized the EEMU.
Finally we made compared CDS prices of each maturity to our TSDPs by term series
regression and found that CDS prices were well explained by DP levels and slopes of
TSDPs. Since the CDS prices are formed in a different market by different players, this
result will show that our approach and model to deriving TSDP measures are effective. In
addition the regression model will enable us to use for trading CDSs.
Overall, our empirical model analysis on credit risk of the main states of the EEMU will
be effective and the results therein will be useful for decision making in credit investment
and risk management.
Optimal Investment in Correlated
Stocks and an Index Bond
for Defined Contribution Pension
2) The price process of index bond with real return r satisfies
dI(t)
= (r + i) dt + σ00 dW0 (t)
I(t)
3) The price process of n stocks are given by
dS (t)
= µdt + Σ dWT (t)
S (t)
Ryou YAMAGATA
Graduate Program of Science and Technology
Mathematics Division
Faculty of Science and Technology, Sophia University
1
where µ is the expected of return vector on the stocks,
Σ is volatility matrix of stocks, and WT (t) is Brownian motion vector with Wk ! W j (k " j, k, j =
0, 1, 2, . . . , n)
Then we denote the market price of risk by θ, which is
described as a function of r, i, µ, R and Σ.
Moreover, define the salary process of the pension plan
member by
Introduction
The purpose of this research work is to derive the optimal investment plan to hedge risks and maximize the total
wealth for defined contribution pension (DCP). The target
of investments are domestic and foreign stocks, domestic
and foreign bonds, Real Estate Investment Trust (REIT),
insurance, trust, fund, and other financial products. We
consider the problem where pension members invest in financial products including stocks and index bonds. On the
assumption that pension members are not day-trader and
their portfolios are rebalanced once a week at most, we derive a theoretical result of the optimal investment plan with
numerical examples.
2
dY(t)
Y(t)
κ = γ + i, σ s
Throughout the following discussion, we assume that
the market is arbitrage-free and complete. Let t ∈ [0, T ]
be continuous time variable provided that t = 0 is the starting time of DC plan, and t = T is the terminal time of the
pension. Then we define the stochastic price level as
=
idt + σ00 dW0 (t)
=
p >0
= σy + σ00 ,
Y(0) = y > 0
π0 (t)
π (t)
1
π(t) = .
..
πn (t)
Basic definition for optimal investment
problem
dP(t)
P(t)
P(0)
κdt + σ s dW(t)
where γ is the expected growth rate of salary, σy is the
volatility of salary, and both are constant. The salary process Y(t) can be solved in the closed form by using Ito’s
Lemma.
Suppose that a DC pension member contributes to an
index bond and n kinds of financial products such as stocks
or funds with the fixed contribute rate c (> 0). Then the
portfolio weight process π(t) is assumed to be
Theory of the optimal investment problem
2.1
=
The first element π0 (t) corresponds to the index bond,
and (π1 (t), π2 (t), . . . , πn (t)) the n kinds of financial prod'
ucts. The member invests share of 1 − ni=0 πi in the riskfree bond. Hence we define the wealth process X π (t) with
an initial value of x (0 ≤ x < ∞) in the following manner.
(
)
dX π (t) = X π (t) Rdt + π% (t)σ (θdt + dW(t)) +cY(t)dt (1)
where cY(t) is the amount of money contributed to the index bond and the financial products at time t. In addition,
we set the stochastic discount factor H(t) by
where the constant i is the expected rate of inflation and
W0 (t) is a Brownian Motion. The volatility of the price
level is σ00 . And we have the following properties:.
*
+
1 2
%
H(t) = exp −Rt − &θ& t − θ W(t)
2
1) The price of risk-free an bond with risk-free rate R is
dB(t)
= Rdt
B(t)
which adjusts for the nominal interest rate and the market
price of risk. The following definition is quoted from [2]
Definition 1 A portfolio vector π(t) is said to be admissible if the corresponding wealth process X π (t) in (1) satisfies
,- T
.
H(s)
π
X (t) + Et
cY(s)ds ≥ 0, P − a.s.
H(t)
t
Then we obtain the total wealth process
where {F (t)}t is the Brownian filtration and Et denotes the
conditional expectation given F (t) under P. We call the
term
,- T
.
H(s)
Et
cY(s)ds
H(t)
t
We consider also that the differential discounted process
of the total wealth. By Ito product rule
Assume that the utility function of the constant relative
risk aversion has the form
2.3
V(t) = X π (t) + D(t)
It follows from (5), (7) and (8) that
dV(t) = dX π (t) + dD(t)
d(H(t)V(t))
=
the human capital.
z1−γ
u(z) =
1−γ
/
0
max E u(X π (T )
subject to
dX π (t) =
(2)
=
π∈Ā(x)
(3)
V ∗ (t) =
(
)
X π (t) Rdt + π% (t)σ (θdt + dW(t))
(4)
x
(5)
with
(9)
x+d
Z1 (t)
H(t)
7
82
1 1−γ
1 − γ %
2
Z1 (t) = exp
θ W(t) −
&θ& t
γ
2
γ
Therefore solving for π∗ (t) gives
8
< =−1 7 1
1 < % =−1
D(t)
∗
%
π (t) =
σ
θ+ σ
θ − σY ∗
γ
γ
V (t) − D(t)
which is the optimal weight process of our aim.
Modeling for optimal investment problem
3
Numerical work
Some numerical works using real financial data have been
presented in the symposium.
We consider the plan member’s future contribution. For
that purpose, define the present value of expected future
contribution process by
,- T
.
H(s)
D(t) = Et
cY(s)ds
H(t)
t
2
c1
=
exp {β(T − t)} − 1 Y(t)
(6)
β
References
[1] Steven E. Shreve. (2004) Stochastic calculus for finance ii continuous-time models, Springer
for t ∈ [0, T ] with
[2] Aihua Zhang. Cristian-Oliver Ewald. (2010) Optimal
investment for a pension fund under inflation risk,
Math Meth Oper Res 71:353-369
β = κ − R − σ s θ0
The representation (6) yields the stochastic differential
equation such as
(
)
dD(t) = D(t) (R + σ s θ0 ) dt + σ%Y dW(t) − cY(t)dt
max E [u(V(T ))]
π∈A1 (x+d)
Thus the problem (3)-(5) is replaced by (9).
Let the superscribe (∗ ) be the corresponding optimal process hereafter, then the optimal wealth process is given by
The problem stated in (3)-(5) is the classical terminal
wealth optimization problem. Our final aim is to derive
the optimal portfolio weight process π(t) as the solution to
the optimization problem.
2.2
Theoretical optimal investment solution
/
0
max E u(X π (T )) =
+cY(t)dt
X π (0)
H(t)dV(t) + V(t)dH(t) + dV(t)dH(t)
In view of the observation in the previous subsection, we
restate the classical terminal wealth optimization problem
(3)-(5) in terms of V(t).
Then, as is discussed in [3], the representative pension plan
member wishes to maximize the expected utility from the
terminal value of the pension fund given an initial investment of x > 0, that is
π∈Ā(x)
(8)
[3] Nicole El Karoui, Monique Jeanblanc-Picque. (1998)
Optimization of consumption with labor income, Finance Stochastic 2, 409-440
(7)
1
•
•
•
–
1
–
CB
5
CB
್ᖜไ㝈䛾⤌䜏
䝃䞊䜻䝑䝖䝤䝺䞊䜹䞊ไᗘ䛾⤌䜏
CBⓎື
ᣑ್ᖜ
ᖒධ್ẁ
❧㛤
2
2.1
∆Xt
∆Xt = φ∆Xt−1 + εt ,
εt
i.i.d.
∼
N (0, σ 2 ), t = 1, 2, . . . .
• X 0 ≡ x0
• σ > 0, |φ| < 1
CB
→ ∆Xt
∆Xt
2
(φ, σ )
→
2.2
•
n
{Pt }nt=1
X 0 ≡ x0
σ2
•
σ 2 (0), φ(0)
φ
• {∆Xt }
X ∼ f ( · |x0 , σ 2 , φ)
•
σ2
{∆Xt }
σ 2 ∼ π1 ( · |∆X, φ)
•
φ
φ ∼ π2 ( · |∆X, σ 2 )
,
3, 4
2.3
VaR
+M −1
{σ 2 (i)}N
i=N
σ2
!2 = 1
σ
M
φ
N +M
"−1
1
φ# =
M
σ 2 (i),
i=N
+M −1
{φ(i)}N
i=N
N +M
"−1
φ(i)
i=N
Φ
•
•
Xt
c
100 α %
&
P {Xt ≤ c} = Φ c − φ#t−t0 Pt0
#t−t0
c=φ
∆Xt0 +
= VaRα×100%
)
!2
σ
3
'
()
t−t
0 −1
"
k=0
!2
σ
#2k
φ
t−t
0 −1
"
k=0
*1/2
φ#2k
*1/2
Φ−1 (α)
VaR
ィ⟬⤖ᯝ
1. 䝟䝷䝯䞊䝍䛾᥎ᐃ್
䝟䝷䝯䞊䝍
᥎ᐃ್
Geweke 䛾ᚋศᕸ᮰ุᐃ
ߪଶ
ͳʹǤͷͶ
ͲǤͶʹͺ ൏ ͳǤͶͷ
ȭ
െͲǤͳͲͲ
᮰
ͲǤͲͲͶሺ൏ ͳǤͶͷሻ ᮰
2. VaR
᥎ᐃ䛧䛯䝟䝷䝯䞊䝍䛷1000ᮏ䛾 path䜢Ⓨ⏕
Ⴀᴗ᪥
VaR1%
1㐌ඛ
ʹʹ͵Ǥͻͷ
2㐌ඛ
ͳͷǤʹ
3㐌ඛ
ͳͳͳ͵Ǥʹ
4㐌ඛ
͵͵ͷǤͷͲ
Statistical Arbitrage on Volatility
1
Ngai Hang Chan
Department of Statistics
Chinese University of Hong Kong
Shatin, NT, Hong Kong
[email protected]
Abstract
This talk discusses recent developments of statistical arbitrage and cointegration. By virtue of some of the asymptotic results about co-integration tests,
a pair-trading strategy is constructed from which statistical arbitrage can
be developed. In particular, a statistical arbitrage strategy is constructed
for volatility trading. The talk concludes with some of the applications to
financial data.
1
Workshop in Statistical Applications and Time Series. November 28–30, 2013,
Kanazawa, Japan. Joint work with Mini Pun and Philip Lee. Research supported in
part by grants from HKSAR-RGC-GRF.
" Unobserved Heterogeneity in Longitudinal Data: An Empirical Bayes
Perspective"
By Roger Koenker, University of Illinois at Urbana-Champaign
Abstract
Empirical Bayes methods for Gaussian and binomial compound decision
problems involving longitudinal data are considered. A new convex
optimization formulation of the nonparametric (Kiefer-Wolfowitz) maximum
likelihood estimator for mixture models can be used to construct
nonparametric Bayes rules for compound decisions. The methods are
illustrated with some simulation examples as well as an application to
predicting baseball batting averages. Comparisons with nonparametric
Bayesian methods based on Dirichlet process priors are also provided.
確率優越性の評価基準で平均と分散の双方に順序制約がある場合の２つの正規母
集団の平均の推定とその応用
目白大学・社会学部 張 元宗 慶應大学・理工学部 篠崎 信雄
1．はじめに
制約された母数空間における推定問題が数多く考えられるが、母平均に関する代表的な線形制約条件は次のような
ものが挙げられる。(1) 非負性 (2) 順序制約 (simple order) (3)simple tree order (4) 傘型順序制約である。例えば、
順序制約は、年齢とともに平均値が大きくなると考えられる量（児童の身長など）、薬品の投与量とともに平均的に
大きくなると考えられる反応量などの場合に考えられている。このような場合の統計的推測について、古くから様々
の研究が進められてきているが、1988 年以前の研究については、Barlow et al.(1972) や Robertson et al. (1988)
で詳しく解説されている。その後の発展、特に、点推定及び区間推定については、Silvapulle & Sen(2005) や van
Eeden(2006) のモノグラフィによって解説されている。また、国内では張 (Chang) と篠崎 (Shinozaki) は制約条件を
考慮する最尤推定量（RMLE）による改良問題について研究を推進し、Kubokawa, Tsukuma, Marchand, Perron,
Strawderman らは一般ベイズ推定量の許容性およびミニマックス性などについて精力的に研究を進めている。
ここでは、母平均と分散の双方に順序制約条件がある場合の２つの正規母集団における平均の推定問題を、確率優越
性 (stochastic domination) あるいは平均 2 乗誤差 (MSE) を評価基準として考える。Xij , i = 1, 2, j = 1, . . . , ni は
!ni
平均 µi 、分散 σi2 の正規分布からの独立な観測値とする。ここで、µi と σi2 は共に未知で、X̄i = j=1
Xij /ni , s2i =
!ni
2
2
j=1 (Xij − X̄i ) /(ni − 1) は µi と σi の不偏推定量である。Oono & Shinozaki(2005) が µi の打ち切り型推定量
µ̂OS
= min{X̄1 , µ̂GD },
1
µ̂OS
= max{X̄2 , µ̂GD }
2
(1.1)
を提案し、µ̂OS
が制約条件を無視した µi の最尤推定量 X̄i を MSE で改良するための必要十分条件を与えた。ここ
i
で、µ̂GD は Graybill & Deal(1959) が提案した２つの正規分布の共通母平均の推定量
µ̂GD =
n1 s22
n2 s21
X̄
+
X̄2
1
n1 s22 + n2 s21
n1 s22 + n2 s21
である。しかし、分散にも順序制約条件 σ12 ≤ σ22 がある場合、s21 > s22 のとき、n1 s22 /(n1 s22 + n2 s21 ) < n1 /(n1 + n2 )
となり不自然である。そこで、次のような推定量
µ̂CS
1
µ̂CS
2
=
=
"
"
µ̂OS
1 ,
1
min{X̄1 , n1n+n
X̄1 +
2
µ̂OS
2 ,
1
max{X̄2 , n1n+n
X̄1 +
2
n2
n1 +n2 X̄2 },
if s21 ≤ s22
if s21 > s22 ,
(1.2)
if s21 ≤ s22
if s21 > s22
(1.3)
n2
n1 +n2 X̄2 },
を考え、µ̂CS
と µ̂OS
との比較を行い、次のような結果を得た。
i
i
２．結果
大きい分散に対応する母平均 µ2 の推定：
定理１． µ̂CS
は µ̂OS
より確率的に優れている。つまり、すべての µ1 ≤ µ2 , σ12 ≤ σ22 に対して、
2
2
"
#
"
#
OS
P | µ̂CS
−
µ
|≤
d
≥
P
|
µ̂
−
µ
|≤
d
,
2
2
2
2
∀d>0
が成立。（証明を付録に置く。）
系２． MSE を基準にしたとき、µ̂CS
は µ̂OS
より優れている。
2
2
小さい分散に対応する母平均 µ1 の推定：
定理３． ∆ = µ2 − µ1 が十分大きいとき、µ̂CS
は µ̂OS
よりも大きな MSE をもつ。
1
1
次に、(µ1 , µ2 ) の同時推定を考える。
CS
OS
OS
2
2
定理４． (µ̂CS
1 , µ̂2 ) は (µ̂1 , µ̂2 ) より次の意味で優れている、つまり、 すべての µ1 ≤ µ2 , σ1 ≤ σ2 に対して、
P
"$
2 %
i=1
µ̂i − µi
σi2 /ni
&2
≤d
#
≥P
"$
2 %
i=1
µ̂OS
− µi
i
σi2 /ni
&2
#
≤d ,
∀ d > 0,
が成立。
参考文献：(1)Chang Y.-T. and Shinozaki, N. (2012) “Estimation of Ordered Means of Two Normal Distributions with Ordered
Variances”, Journal of Mathematics and System Science Vol.2, No.1, pp. 1-7. (2)Chang Y,-T. Oono, Y. and Shinozaki, N. (2012),
“Improved estimators for the common mean and ordered means of two normal distributions with ordered variances”, Journal of
Statistical Planning and Inference, 142, 2619-2628. (3)Graybill, F.A. and Deal, R. B. (1959), “Combining unbiased estimators,”
Biometrics 15, pp. 543-550. (4)Hwang, J. T. (1985), “Universal domination and stochastic domination,” Ann. Statist., Vol.13,
No.1, pp. 295-314. (5)Oono, Y. and Shinozaki, N. (2005), “Estimation of two order restricted normal means with unknown and
possibly unequal variances,” Journal of Statistical Planning and Inference, Vol. 131, Issue2, pp. 349-363. (6) Silvapulle, M. J.,
Sen, P. K. (2004). Constrained Statistical Inference. Wiley, New Jersey. (7) van Eeden, C., (2006). Restricted Parameter Space
Estimation Problems. Lecture notes in Statistics 188, Springer.
謝辞：本研究は科研費 (基盤研究 (C)No.22500263) 助成を受けたものである。また、著者張元宗は目白大学の特別研究費の助成も受けています。
付録 （定理１の証明）
まず、P {|µ̂CS
2 − µ2 | ≤ d|X̄1 > X̄2 } を下記のように表現できる。
"
#
"
#
2
2
2
2
CS
2
2
2
2
P |µ̂CS
−
µ
|
≤
d|
X̄
>
X̄
,
s
≤
s
P
{s
≤
s
}
+
P
|µ̂
−
µ
|
≤
d|
X̄
>
X̄
,
s
>
s
2
1
2 1
2
1
2 1
2
2
1
2
2
2 P {s1 > s2 }.
µ̂CS
は µ̂OS
より、確率的に優れているために、次の式を証明すればよい。
2
2
"
#
"
#
2
2
OS
2
2
P |µ̂CS
−
µ
|
≤
d|
X̄
>
X̄
,
s
>
s
≥
P
|µ̂
−
µ
|
≤
d|
X̄
>
X̄
,
s
>
s
2
1
2 1
2
1
2 1
2
2
2
2 .
c = n1 /(n1 + n2 ) とする。 (A.1) の左辺は
"
#
2
2
P |c(X̄1 − µ2 ) + (1 − c)(X̄2 − µ2 )| ≤ d|X̄1 > X̄2 , s1 > s2 .
(A.1)
(A.2)
になる。次に、変数変換 Wi = X̄i − µ2 , i = 1, 2. を行うと、W1 ∼ N (−∆, τ12 ), W2 ∼ N (0, τ22 ) に従い、W1 と W2
は互いに独立になる。ここで、∆ = µ2 − µ1 であり、τi2 = σi2 /ni , i=1,2 である。
よって、(A.2) は
P {−d ≤ cW1 + (1 − c)W2 ≤ d|W1 > W2 , s21 > s22 }.
(A.3)
になる。さらに、変数変換 V1 = W1 − W2 , V2 = (τ22 /τ12 )W1 + W2 . を行うと、V1 と V2 は互いに独立で、V1 ∼
N (−∆, τ12 + τ22 ), V2 ∼ N (−(τ22 /τ12 )∆, τ22 (τ12 + τ22 )/τ12 ) である。また、W1 =
ら、便宜上、s21
>
s22
を固定し、(A.3) は
τ12 (V1 +V2 )
, W2
τ12 +τ22
=
τ12 V2 −τ22 V1
,
τ12 +τ22
であるか
"
#
"
#
(cτ12 − (1 − c)τ22 )V1 + τ12 V2
P −d ≤ cW1 + (1 − c)W2 ≤ d|W1 > W2 = P −d ≤
≤
d|V
>
0
1
τ12 + τ22
"
#
"
#
= P (cτ12 − (1 − c)τ22 )V1 + τ12 V2 ≥ −d(τ12 + τ22 )|V1 > 0 − P (cτ12 − (1 − c)τ22 )V1 + τ12 V2 ≥ d(τ12 + τ22 )|V1 > 0 .
= E{g(c, V1 )|V1 > 0},
になる。ここで、
%
((1 − c)τ22 − cτ12 )v1 + τ22 ∆ + d(τ12 + τ22 )
'
g(c, v1 ) = Φ
τ1 τ2 τ12 + τ22
&
%
&
((1 − c)τ22 − cτ12 )v1 + τ22 ∆ − d(τ12 + τ22 )
'
−Φ
.
τ1 τ2 τ12 + τ22
である。同様に、(A.1) の右辺については、次のように計算される。
"
#
OS
2
2
P |µ̂2 − µ2 | ≤ d|X̄1 > X̄2 , s1 > s2 = E{g(γ0 , V1 )|V1 > 0, s21 > s22 },
ここで、γ0 = n1 s22 /(n1 s22 + n2 s21 ) である。
v1 > 0 を固定した上で、g(u, v1 ) は 0 ≤ u ≤ 1 の非減少関数であることを簡単に証明できる。また、s21 > s22 に対
して、c > γ0 が成立するので、
"
#
"
#
2
2
P |µ̂2 − µ2 | ≤ d|X̄1 > X̄2 , s21 > s22 ≥ P |µ̂OS
−
µ
|
≤
d|
X̄
>
X̄
,
s
>
s
2
1
2 1
2
2 ,
が成立、証明が完成した。
LMSR-method
M2
1
LMS-method
Cole and Green (1992)
(1990),Cole,Freeman and Preece(1998)
) LMS-mehod
3
1
LMS-method
LMSR-method
LMSR-method
1
2
2.1
(Cole
LMSR-method
LMS-method
{(xi , yi ) : i = 1, · · · , n}
zi =
xi
yi
L(xi ), M (xi ), S(xi )(i = 1, · · · , n)
%
&L(xi )
yi
−1
M (xi )
L(xi )S(xi )
, L(xi ) "= 0
%
&
yi
log
M (xi )
S(xi )
3
, i = 1, · · · , n
, L(xi ) = 0
i.i.d.
Z1 , · · · , Zn ∼ N (0, 1)
L = (L(x1 ), · · · , L(xn ))T , M = (M (x1 ),· · · , M (xn ))T ,
S = (S(x1 ), · · · , S(xn ))T
L, M , S Natural Cubic Spline(NCS)
Π(L, M , S)
αL T
αM T
αS T
Π(L, M , S) = !(L, M , S) −
L KL −
M KM −
S KS
2
2
2
αL , α M , αS
!(L, M , S)
!(L, M , S) =
n (
'
L(xi ) log
i=1
1
%
yi
M (xi )
&
− log S(xi )
)
K
xi (i = 1, · · · , n)
n×n
(Green and Silverman(1994)
)
Π(L, M , S)
L, M , S
+
,T
+
,T
+
,T
* = L(x
- = M
* = S(x
* 1 ), · · · , L(x
* n) , M
-(x1 ), · · · , M
-(xn ) , S
* 1 ), · · · , S(x
* n)
L
.+
,
/ .+
,
/ .+
,
/
* i ) : i = 1, · · · , n , xi , M
-(xi ) : i = 1, · · · , n , xi , S(x
* i ) : i = 1, · · · , n
xi , L(x
3
2.2
{(xi , yi ) : i = 1, · · · , n}
NCS
*
-(x), S(x)
*
L(x),
M
LMSR-method
LMS-method
LMSR-method
2
2.2.1
{(xi , y i ) : i = 1, · · · , n}, y i = (yi1 , · · · , yim )T
Lk (xi ), Mk (xi ), Sk (xi )
xi
yik
3
zik =
%
&Lk (xi )
yik
−1
Mk (xi )
Lk (xi )Sk (xi )
%
&
y
ik
log
Mk (xi )
Sk (xi )
, Lk (xi ) "= 0
, i = 1, · · · , n; k = 1, · · · , m
, Lk (xi ) = 0
Z i = (Zi1 , · · · , Zim )T
Z i ∼ Nm (0, I(R(xi ))), Z i ⊥
⊥ Z j (i "= j , i, j = 1, · · · , n)
I(R(xi ))
R(xi )
T
Lk = (Lk (x1 ), · · · , Lk (xn )) , M k = (Mk (x1 ), · · · , Mk (xn ))T , S k = (Sk (x1 ), · · · , Sk (xn ))T
(k = 1, · · · , m)
R = (R(x1 ), · · · , R(xn ))T
Lk , M k , S k (k = 1, · · · , m)
R NCS
T
T T
T
T T
T
L = (L1 , · · · , Lm ) , M = (M 1 , · · · , M m ) , S = (S 1 , · · · , S Tm )T
Π(L, M , S, R)
Π(L, M , S, R) = !(L, M , S, R) −
αL T
αM T
αS T
αR T
L KL −
M KM −
S KS −
R KR
2
2
2
2
2
αL , α M , α S , α R
!(L, M , S, R)
0m (
%
n
'
'
!(L, M , S, R) =
Lk (xi ) log
i=1
k=1
yik
Mk (xi )
&
− log Sk (xi )
1
1
− log |I(R(xi ))| − z Ti I(R(xi ))−1 z i
2
2
)
)
z i = (zi1 , · · · , zim )T (i = 1, · · · , n)
K
nm × nm
K =
Im ⊗ K
Im m × m
⊗
Π(L, M , S, R)
Lk , M k , S k (k =
T
*
*
*
-k (xn ))T
1, · · · , m), R
Lk = (Lk (x1 ), · · · , Lk (xn ))
M k = (Mk (x1 ), · · · , M
* k = (S*k (x1 ), · · · , S*k (xn ))T (k = 1, · · · , m), R
* = (R(x
* 1 ), · · · , R(x
* n ))T
S
.+
,
/ .+
,
/ .+
,
/
* k (xi ) : i = 1, · · · , n , xi , M
-k (xi ) : i = 1, · · · , n , xi , S*k (xi ) : i = 1, · · · , n ,
xi , L
k = 1, · · · , m,
.+
,
/
* i ) : i = 1, · · · , n
xi , R(x
{(xi , y i ) : i = 1, · · · , n}, y i = (yi1 , · · · , yim )T
* k (x), M
-k (x), S*k (x)(k = 1, · · · , m), R(x)
*
L
NCS
3m + 1
2.2.2
i ∈ {1, · · · , n}
R(xi ) = ρ
i.i.d.
Z 1 , · · · , Z n ∼ Nm (0, I(ρ))
I(ρ)
Π(L, M , S, ρ)
ρ
Π(L, M , S, ρ) = !(L, M , S, ρ) −
αL T
αM T
αS T
L KL −
M KM −
S KS
2
2
2
αL , α M , αS
!(L, M , S, ρ) =
!(L, M , S, ρ)
0m (
n
'
'
%
&
)
yik
Lk (xi ) log
− log Sk (xi )
Mk (xi )
i=1 k=1
)
1
1
− log |I(ρ)| − z Ti I(ρ)−1 z i
2
2
z i = (zi1 , · · · , zim )T (i = 1, · · · , n)
Π(L, M , S, ρ)
T
* k = (L
- k = (M
* k (x1 ), · · · , L
* k (xn )) , M
-k (x1 ), · · · , M
-k (xn ))T ,
L
* k = (S*k (x1 ), · · · , S*k (xn ))T (k = 1, · · · , m)
S
.+
,
/ .+
,
/ .+
,
/
* k (xi ) : i = 1, · · · , n , xi , M
-k (xi ) : i = 1, · · · , n , xi , S*k (xi ) : i = 1, · · · , n ,
xi , L
3
k = 1, · · · , m
NCS
3m
{(xi , y i ) : i = 1, · · · , n}, y i = (yi1 , · · · , yim )T
* k (x), M
-k (x), S*k (x)(k = 1, · · · , m) ρ*
L
2.2.3
NCS
-k (x)(k = 1, · · · , m)
M
- k = (M
-k (x1 ), · · · , M
-k (xn ))T (k = 1, · · · , m)
M
.+
,
/
-k (xi ) : i = 1, · · · , n , k = 1, · · · , m
xi , M
Berwin(2005,Sec.2.2)
NCS
NCS
-k (x)(k = 1, · · · , m)
M
*
R(x)
*
−1 ≤ R(x)
≤1
Berwin(2005,Sec.2.3)
* = (R(x
* 1 ), · · · , R(x
* n ))T
R
.+
,
/
* i ) : i = 1, · · · , n
xi , R(x
NCS
*
R(x)
[1] Cole,T.J. and Green,P.J.(1992). Smoothing reference centile curves: the LMS method
and penalized likelihood,Statistics in Medicine,11,1305-1319.
[2] Cole,T.J.(1990). The LMS method for constructing normalized growth standards,
European Journal of Clinical Nutrition,44,45-60.
[3] Cole,T.J.,Freeman,J.V. and Preece,M.A.(1998). British 1990 growth reference centiles for weight, height, body mass index and head circumference fitted by maximum
penalized likelihood.Statistics in Medicine,17,407-429.
[4] Green,P.J. and Silverman,B.W.(1994). Nonparametric regression and generalized linear models: a roughness penalty approach,Vol.58 of Monographs on Statistics and
Applied Probability,Chapman and Hall,London.
[5] Berwin,A.T.(2005). Shape constrained smoothing using smoothing splines,Computational Statistics,20,81-103.
4
ศᩓඹศᩓ⾜ิᇶ࡙ࡃุูศᯒࡢ⤊↉
᪂ᮧ ⚽୍͊
͊ᡂ㋣Ꮫ⤒῭Ꮫ㒊ࠛ ᮾி㒔Ṋⶶ㔝ᕷྜྷ⚈ᑎ⏫ 㸬 ྜ ྰ ุ ᐃ ࡼ ࡿ 㔜 ᕪ ᳨ ド ἲ ࡣ ࡌ ࡵ ホ ౯ ⏝ ࠸ ࡓ ࢹ ࣮ ࢱ ᭱ ᑠ ㄗ ศ 㢮 ᩘ 00 10 ᇶ ࡙ ࡃ ᭱ 㐺 ⥺ ᙧ ุ ู
010 ࡢ ࢹ ࣮ ࢱ ࡋ ࡚ ᖺ ࡢ ࠕ ⤫ ィ ධ 㛛 ࠖ
2 /')ࢆ 㛤 Ⓨ ࡋ ࡓ 㸬 010 ᇶ ‽ ᇶ ࡙ ࡃ ⥺ ᙧ
㛵 ᩘ 2
ࡢ ୰ 㛫 ヨ 㦂 ࡢ ࢹ ࣮ ࢱ ࢆ ⏝ ࠸ ࡿ >@㸬 ᢥ ၥ
ุู㛵ᩘࡣΰྜᩚᩘィ⏬ἲ࡛ࡋᐇ⌧࡛ࡁࡎ㸪
ࡢ ヨ 㦂 ࢆ 7 ࡽ 7 ࡢ ၥ ၥ ศ ࡅ ࡓ 㸬ࡇ ࡢ
᭱ ึ ࡣ ,32/') ࡪ ࣔ ࢹ ࣝ ࢆ 㛤 Ⓨ ࡋ 㸪)LVKHU
ಶ ࡢ ᚓ Ⅼ ࢆ ㄝ ᫂ ኚ ᩘ ࡋ 㸪ྜ ィ ᚓ Ⅼ ࡢ 㸣 Ⅼ
ࡢ ࣜ ࢫ ࢹ ࣮ ࢱ 㸪 ࢫ ࢫ 㖟 ⾜ ࢹ ࣮ ࢱ 㸪 &3' ࢹ
Ⅼ 㸪 㸣 Ⅼ Ⅼ 㸪 㸣 Ⅼ Ⅼ ࡢ Ỉ
࣮ ࢱ 㸪Ꮫ ⏕ ࢹ ࣮ ࢱ ࢆ ◊ ✲ ࢹ ࣮ ࢱ ࡋ ࡚ ⏝ ࠸ ࡚ ◊
‽ ࡛ ྜ ྰ ุ ᐃ ࢆ ⾜ ࠺ 㸬㸣 Ⅼ ࡛ ࡣ ၥ ၥ 㸪㸣
✲ ࢆ ⾜ ࡗ ࡓ 㸬 ࡑ ࡋ ࡚ 㸪 㸧 ࢫ ࢫ 㖟 ⾜ ⣬ ᖯ ࢹ ࣮
Ⅼ 㸣 Ⅼ ࡛ ࡣ ၥ ၥ ࡛ 010 ࡛ ࠶ ࡿ 㸬 ಶ
ࢱ ࡀ ኚ ᩘ ࡛ ⥺ ᙧ ศ 㞳 ྍ ⬟ ࡛ ࠶ ࡾ 010 ࡢ ༢ ㄪ ῶ
ࡢ ㄝ ᫂ ኚ ᩘ ࡢ ุ ู ࣔ ࢹ ࣝ ࡣ 㒊 ࡛ ಶ ࠶ ࡿ ࡀ 㸪
ᑡ ᛶ ࡽ ࡑ ࡢ ኚ ᩘ ࢆ ྵ ࡴ ಶ ࡢ ุ ู ࣔ ࢹ ࣝ
ኚ ᩘ ࡢ ุ ู ࡣ ព ࡀ ࡞ ࠸ ࡢ ࡛ ኚ ᩘ ௨ ୖ ࡢ ࡀ ⥺ ᙧ ศ 㞳 ྍ ⬟ ࡛ ࠶ ࡿ 㸬ุ ู ಀ ᩘ ㄗ ศ 㢮 ᩘ
ಶ ࡛ ᳨ ウ ࡍ ࡿ 㸬 010 ࡢ ༢ ㄪ ῶ ᑡ ᛶ ࡽ 㸪 㸣 Ⅼ
ࡢ 㛵 ಀ ࡀ ศ ࡾ 㸪ุ ู ಀ ᩘ ࡢ ✵ 㛫 ࡛ ᭱ 㐺 ฝ య ࡢ
ࡢ ྜ ྰ ุ ᐃ ࡛ ࡣ ಶ ࡢ ุ ู ࣔ ࢹ ࣝ ࡀ 010 ࡛ 㸪
ෆ Ⅼ ࢆ ┤ ᥋ ồ ࡵ ࡿ ᨵ ᐃ ,32/') ࢆ 㛤 Ⓨ ࡋ ࡓ 㸬ᮏ
ಶ ࡀ 010 ࡛ ࡞ ࠸ 㸬㸣 Ⅼ 㸣 Ⅼ ࡛ ࡣ ಶ
ᡭ ἲ ࡣ 010 ࡢ ࢹ ࣮ ࢱ ࢆ ⌮ ㄽ ⓗ ุ ู ࡛ ࡁ ࡿ ၏
ࡔ ࡅ ࡀ 010 ࡛ ಶ ࡀ 010 ࡛ ࡞ ࠸ 㸬 ࡇ ࢀ
୍ ࡢ ุ ู ᡭ ἲ ࡛ ࠶ ࡾ 㸪ุ ู ㉸ ᖹ 㠃 ୖ ࡢ ࢣ ࣮ ࢫ ࡢ
ࡼ ࡗ ࡚ 㸪 010 ࢆ ྵ ࡴ ุ ู ศ ᯒ ࡢ ၥ 㢟 ࡀ ᳨ ウ ࡛
ุ ู ࢆ ṇ ࡋ ࡃ ฎ ⌮ ࡀ ࡛ ࡁ ࡿ 㸬ୖ グ ࡢ ✀ ࡢ ᐇ
ࡁ ࡿ 㸬 690 ࡣ +690 ࡛ 010 ࡞ ࢹ ࣮ ࢱ ࡢ ุ ู ࢆ
ࢹ ࣮ ࢱ ࡽ ࢣ ࣮ ࢫ ᩘ ಸ ࡢ ࣜ ࢧ ࣥ ࣉ ࣜ ࣥ ࢢ ᶆ
ุูศᯒࡢฟⓎⅬࡋࡓࡇࡣ㠉᪂ⓗ࡛࠶ࡿ㸬
ᮏ ࢆ స ᡂ ࡋ 㸪 010 ࡢ ㏆ ఝ ್ ࢆ 㧗 ㏿ ࡛ ồ ࡵ ࡿ ᨵ ᐃ
ࡋ ࡋ 㸪 690 ࢆ ྵ ࡵ ࡇ ࢀ ࡲ ࡛ ุ ู ศ ᯒ ࡣ 010 ,3/32/') ࢆ /') ࣟ ࢪ ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ ࡛ ࡢ ࢹ ࣮ ࢱ ࡢ ᐇ ド ◊ ✲ ࢆ ⾜ ࡗ ࡚ ࡇ ࡞ ࡗ ࡓ 㸬ุ ู
㔜 ᕪ ᳨ ド ἲ ࢆ ⏝ ࠸ ࡚ ẚ ㍑ ࡋ 㸪 ᨵ ᐃ ,3/32/')
ศᯒࡢ㔜せ࡞ᛂ⏝ศ㔝࡛࠶ࡿࣃࢱ࣮ࣥㄆ㆑࡛
ࡀᅽಽⓗᏛ⩦ᶆᮏ᳨ドᶆᮏ࡛Ⰻ࠸⤖ᯝࢆ
ࡣ 㸪 010 ࡢ ࢹ ࣮ ࢱ ࡢ ุ ู ࡀ 㔜 せ ࡛ ࠶ ࡿ 㸬 ࡇ ࡢ
ᚓ ࡓ 㸬 ௨ ୖ ࡢ ⤖ ᯝ ࢆ ㋃ ࡲ ࠼ ᖺ ◊ ✲ ᡂ ᯝ
ࡓ ࡵ 㸪 ᅇ 010 ࡀ ⮬ ᫂ ࡞ ヨ 㦂 ࡢ ྜ ྰ ุ ᐃ ࢹ ࣮
ࢆ ࡲ ࡵ ࡚ ฟ ∧ ࡋ ࡓ 㹙 㹛㸬 ࢱ ࢆ ⏝ ࠸ ࡿ ࡇ ࡣ 㐺 ษ ⪃ ࠼ ࡿ 㸬ヨ 㦂 ࡢ ྜ ྰ ุ
ᖺ ࡽ ᛂ ⏝ ◊ ✲ ࡋ ࡚ 010 ࡢ ุ ู ↔
ᐃ ࢹ ࣮ ࢱ ࡣ 㸪ከ ࡃ ࡢ ே ࡀ ධ ᡭ ࡋ ᐜ ᫆ ᳨ ド ࡛ ࡁ 㸪
Ⅼࢆ⤠ࡗ࡚◊✲ࢆ⾜࠺ࡇࡋࡓ㸬ࡋࡋ
ྜ᱁Ỉ‽ࢆ✀ࠎኚ࠼ࡿࡇ࡛ ⩌ࡢࢣ࣮ࢫᩘ
010 ࡢ ࢹ ࣮ ࢱ ࡣ ⌧ ᐇ ከ ࡃ ࡞ ࠸ ࡀ 㸪 タ ၥ ࡢ ᚓ
ࡢ✀ࠎࡢྜࡢ㐪࠸ࡼࡿุูࡢᙳ㡪ࢆ᳨
Ⅼࢆㄝ᫂ኚᩘࡋ࡚ᚓⅬྜィ࡛ྜྰุᐃࢆ⾜
ウ࡛ࡁࡿ࡞ࡢⅬࡀ࠶ࡿ㸬
࠺ ⩌ ุ ู ࡣ 㸪⮬ ᫂ ࡞ 010 ࡞ ࡿ ⥺ ᙧ ุ ู 㛵
㸣 Ỉ ‽ ࡢ 㔜 ᕪ ᳨ ド ἲ ᩘ ࡀ ồ ࡲ ࡿ 㸬 ಶ ࡢ ྜ ྰ ุ ᐃ ࡛ /') ࡢ ㄗ ศ 㢮 ☜
⾲ ࡣ 㔜 ᕪ ᳨ ド ἲ ࡢ ⤖ ᯝ ࡛ ࠶ ࡿ 㸬 ิ
⋡ ࡢ ⠊ ᅖ ࡣ >@㸪4') ࡣ >@
┠ࡣྛุูᡭἲࡢୗࣔࢹࣝ␒ྕࢆ♧ࡍ㸬 ิ
ࡁ ࠸ ࡇ ࡀ ศ ࡗ ࡓ 㸬 ࡲ ࡓ 㸪 4') ࡣ ୍ ᪉ ࡢ ⩌ ࡢ
┠ ࡣ ㄝ ᫂ ኚ ᩘ ࡢ ᩘ ࢆ ⾲ ࡍ 㸬 ิ ࡣ Ꮫ ⩦ ᶆ ᮏ ࡢ ኚ ᩘ ್ ࡀ ୍ ᐃ ್ ࡢ ሙ ྜ ࡢ ࡳ 㸪୍ ⯡ ㏫ ⾜ ิ ࢆ ⏝
ಶ ࡢ ᖹ ᆒ ㄗ ศ 㢮 ☜ ⋡㸦 0($1㸧࡛ ࠶ ࡿ 㸬 ิ ࡣ ᳨
࠸ ࡓ -03㹙 㹛 ࡢ 4') ࡸ ṇ ๎ ุ ู ศ ᯒ ࡀ ⩌
ド ᶆ ᮏ ࡢ ᖹ ᆒ ㄗ ศ 㢮 ☜ ⋡ 㸦 0($1㸧 ࡛ 㸪 ᭱ ᑡ ࡢ
ᣦᐃࡋࡓࢣ࣮ࢫࢆ ⩌ㄗุูࡍࡿ㔜࡞
ࣔ ࢹ ࣝ ࢆ 㔜 ᕪ ᳨ ド ἲ ࡀ 㑅 ࡪ ࣔ ࢹ ࣝ ⪃ ࠼
⍗ ⑅ ࡢ ཎ ᅉ ࡀ ศ ࡗ ࡓ >@㸬 ࡲ ࡓ ࠕ ᑠ ᶆ ᮏ ࡢ
ࡿ㸬
ࡓ ࡵ ࡢ N㔜 ᕪ ᳨ ド ἲ ࠖࢆ 㛤 Ⓨ ࡋ 㸪ୖ グ ✀ ࡢ
ᨵ ᐃ ,32/')㸦 ⾲ ࡢࠕ ,3ࠖ㸧ࡣ 㔜 ᕪ ᳨ ド
ࢹ ࣮ ࢱ ຍ ࠼ ࡚ 㸪 ᖺ ᗘ ࡢ ➹ ⪅ ࡢ ⤫ ィ ධ 㛛 ࡢ
ἲࡢᏛ⩦ᶆᮏ࡛ඖࢹ࣮ࢱྠࡌࡃ ಶࡢุูࣔ
୰ 㛫 ヨ 㦂 ࡛ 㔜 ᕪ ᳨ ド ἲ ࢆ ⾜ ࠸ 㸪 ᨵ ᐃ ,3
ࢹ ࣝ ࡛ 010 ࡛ ࠶ ࡿ ࡀ 㸪 ᳨ ド ᶆ ᮏ ࡛ ࡣ 61 2/') ࡣ Ꮫ ⩦ ᶆ ᮏ ࡛ 㐣 Ꮫ ⩦ ࡍ ࡿ ࡀ 㸪 ᳨ ド ᶆ ᮏ ࡛
ࡔࡅࡀ ࡞ࡗࡓ㸬
ࠕ ᕪ ࠖࡣࠕ 0($1ࠖࠕ 0($1ࠖ
/')ࣟ ࢪ ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ 㸪ࢯ ࣇ ࢺ ࣐ ࣮ ࢪ ࣥ ᭱
ࡢ ᕪ ࡛ ࠶ ࡿ 㸬୍ ⯡ ⓗ 㸪ᨵ ᐃ ,32/') ࡣ Ꮫ ⩦
690㸦 6690㸧 ẚ ࡚ ᳨ ド ᶆ ᮏ ࡛ ࡶ ᖹ ᆒ ㄗ ศ
ࢹ ࣮ ࢱ ࡛ 㐣 Ꮫ ⩦ ࡋ 㸪᳨ ド ᶆ ᮏ ࡛ ุ ู ࡢ ண ⤖ ᯝ
㢮☜⋡ࡀ୍␒ᑠࡉ࠸ࡇࡀศࡗࡓ㸬
ࡣ ᝏ ࡃ ࡞ ࡿ ⪃ ࠼ ࡽ ࢀ ࡓ ࡀ 㸪ᐇ 㝿 ࡣ ␗ ࡞ ࡗ ࡚ ࠸ ⾲ Ỉ ‽ ࡢ ㄗ ศ 㢮 ☜ ⋡ ,3 3 0($1 0($1 ᕪ /RJL 3 0($1 0($1 ᕪ 0($1 ᕪ
690 3 0($1 0($1 ᕪ 0($1 ᕪ
/') 3 0($1 0($1 ᕪ 0($1 ᕪ
ࡓ 㸬ࣟ ࢪ ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ 6690 ࡣ Ꮫ ⩦ ᶆ ᮏ ࡛
ඖ ࢹ ࣮ ࢱ ྠ ࡌ ࡃ ಶ ࡢ ุ ู ࣔ ࢹ ࣝ ࡛ 010 ࡛
࠶ ࡾ 㸪᳨ ド ᶆ ᮏ ࡛ 61 ࢆ 㑅 ࢇ ࡔ ࡀ ㄗ ศ 㢮 ᩘ ࡣ ࡛ ࡞ ࠸ 㸬 ᕪ ࡣ ࡛ ᝏ ࡃ ࡞ ࠸ 㸬 /') ࡣ
Ꮫ⩦ᶆᮏ᳨ドᶆᮏ࡛ᖹᆒㄗศ㢮☜⋡ࡀ ࡢࡶ
ࡢ ࡀ ࡞ ࡃ 61 ࢆ 㑅 ࢇ ࡔ ࡀ 㸪 ᖹ ᆒ ㄗ ศ 㢮 ☜ ⋡ ࡣ
㠀 ᖖ ᝏ ࠸ 㸬 0($1 ᕪ ࡣ 㸪 ࣟ ࢪ
ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ 㸪6690㸪/') ࡢ 0($1 ࡽ ᨵ ᐃ
,32/') ࡢ ࡑ ࢀ ࢆ ᘬ ࠸ ࡓ ࡶ ࡢ ࡛ 㸪 ⠊ ᅖ ࡣ >
@>@>@࡛ ࠶ ࡿ 㸬 ࣟ ࢪ
ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ 6690 ࡣ ᨵ ᐃ ,32/') ࡼ ࡾ
㧗 ࠎ ᝏ ࠸ ࡀ 㸪 /') ࡣ ௨ ୖ ᝏ ࠸ 㸬 Ỉ ‽ 㸣 Ỉ ‽ ࡢ ྜ ྰ ุ ᐃ ᚓ Ⅼ ศ ᕸ ࡢ Ⅼ Ⅼ ࡛ ྜ ྰ ุ ᐃ ࢆ ⾜ ࠺ 㸬ࣟ
ࢪ ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ ࡢ ࠕ 0($1 ᕪ ࠖ ࡢ ⠊ ᅖ ࡣ >
@㸪 690 ࡣ >@㸪 /') ࡣ >@
࡛ ࠶ ࡿ 㸬 /') ࡣ ࡢ ᡭ ἲ ࡀ ᭱ ᑠ ࡞ ࡿ ࣔ ࢹ ࣝ ࡛ ᨵ ᐃ ,32/') ࡼ ࡾ 㸣 ࡶ ᝏ ࠸ 㸬 Ⅼ ࡢ ྜ ྰ
ุ ᐃ ࡛ 㸪 ࣟ ࢪ ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ ࡢ ࠕ 0($1 ᕪ ࠖ ࡢ ⠊
ᅖ ࡣ >@࡛ 㸪 690 ࡣ >@࡛ ࠶ ࡿ ࡢ
ᑐ ࡋ ࡚ 㸪 /') ࡣ >@ ᫂ ࡽ ᝏ ࠸ 㸬
ࣟ ࢪ ࢫ ࢸ ࢵ ࢡ ᅇ ᖐ ࡢ ᖹ ᆒ ㄗ ศ 㢮 ☜ ⋡ ࡣ 㸪ࣔ ࢹ ࣝ
㑅 ᢥ ↓ 㛵 ಀ ࡞ ␒ ┠ ࡢ ࣔ ࢹ ࣝ ࡛ ᨵ ᐃ ,32/') ࡼ
ࡾᑡ࡞࠸㸬
ࡲ ࡵ ᮏ ◊ ✲ ࡛ ࡣ 㸪ヨ 㦂 ࡢ ྜ ྰ ุ ᐃ ࢹ ࣮ ࢱ ࢆ ⏝ ࠸ ࡚
010 ࡢ ุ ู ศ ᯒ ࡢ ၥ 㢟 ࢆ ᳨ ウ ࡋ ࡓ 㸬010 ᇶ ‽
ᇶ ࡙ ࡃ ᨵ ᐃ ,32/') ࡣ Ꮫ ⩦ ᶆ ᮏ ࡛ 㐣 Ꮫ ⩦ ࡋ 㸪
᳨ ド ᶆ ᮏ ࡛ ỗ ⬟ ຊ ࡀ ᝏ ࠸ ⪃ ࠼ ࡽ ࢀ ࡿ ࡀ 㸪ᐇ
㝿 ࡣ ṇ つ ศ ᕸ ࢆ ௬ ᐃ ࡍ ࡿ /') ࡀ Ꮫ ⩦ ᶆ ᮏ ᳨
ド ᶆ ᮏ ࡛ ☜ ⋡ ศ ᕸ ࢆ ௬ ᐃ ࡋ ࡚ ࠸ ࡞ ࠸ 6690 ࡸ ࣟ
ࢪࢫࢸࢵࢡᅇᖐࡼࡾࡶᝏࡗࡓ㸬ࡇࢀࡣ
)LVKHU ࡢ ௬ ㄝ ࢆ ‶ ࡓ ࡍ ⌧ ᐇ ࡢ ࢹ ࣮ ࢱ ࡀ ᑡ ࡞ ࠸
ࡢࡑࢀ࡛⌮ㄽᵓ⠏ࡋࡓࡓࡵ⪃࠼ࡽࢀࡿ㸬
ᩥ ⊩ >@᪂ ᮧ ⚽ ୍ 㸪-03 ࡼ ࡿ ⤫ ィ Ꮫ ࡗ ࡚ ࠾ ࡁ ຮ ᙉ ἲ 㸪 ᮾ ி 㸪 ㅮ ㄯ ♫ 㸪 㸬 >@᪂ ᮧ ⚽ ୍ ([FHO /,1*2 ࡛ Ꮫ ࡪ ᩘ ⌮ ィ ⏬ ἲ ၿ 㸬 >@᪂ ᮧ ⚽ ୍ ᭱ 㐺 ⥺ ᙧ ุ ู 㛵 ᩘ 㸬 ᪥ ⛉ ᢏ 㐃 ฟ ∧ ♫ >@᪂ ᮧ ⚽ ୍ ᩘ ⌮ ィ ⏬ ἲ ࡼ ࡿ ၥ 㢟 ゎ Ỵ ἲ 㸬᪥ ⛉ ᢏ 㐃 ฟ ∧ ♫ >@᪂ ᮧ ⚽ ୍ ⤫ ィ ᩍ ⫱ ࠾ ࡅ ࡿ ุ ู ศ ᯒ ࡢ ᨵ ၿ Ⅼ 㸬
⤫ ィ ᩍ ⫱ ᐇ ㊶ ◊ ✲ ➨ ᕳ ⤫ ィ ᩘ ⌮ ◊ ✲ ᡤ ◊ ✲ ࣞ ࣏ ࣮ ࢺ 㸪 ࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࢆ⏝࠸ࡓ㠀⥺ᙧุูᡭἲࡢᥦ
ᒣᮏࡅ࠸Ꮚ ภ㤋ᕤᴗ㧗➼ᑓ㛛Ꮫᰯ ᐮἙỤ㞞ᙪ 㔠ἑᏛ⤒῭Ꮫ㢮
㸯㸬ࡣࡌࡵ
ࢡࣛࢫูṇゎศᕸࢆసᡂࡋ࡚࠾ࡃࡇ࡛㸪ᐇ㝿
ࣃࢱ࣮ࣥㄆ㆑࠾ࡅࡿ㠀⥺ᙧุูၥ㢟ࡣ㸪ࢧ࣏
ࡢุูࡣ㸪ุูᑐ㇟ࢹ࣮ࢱࡢศᕸ᥎ᐃ㢮ఝ
࣮ࢺ࣋ࢡࢱ࣮࣐ࢩࣥ㸦SVM㸧[1]ࡢฟ⌧ࡼࡾ㸪
㔞ࡢ⟬ฟࡢࡳ࡛ࢡࣛࢫุูࢆ⾜࠺ࡇࡀ࡛ࡁࡿࡓ
ࡁ࡞㐍ᒎࢆࡆ࡚࠸ࡿ㸬ᡃࠎࡣ㸪࣮࢝ࢿࣝᐦᗘ᥎
ࡵ㸪ຠ⋡ⓗࡘỗ⏝ⓗ࡞ᡭἲ࠸࠼ࡿ㸬
ᐃἲ㸦௨㝆㸪KDE ␎ࡍ㸧ࢆ⏝࠸ࡓ㠀⥺ᙧุูᡭ
ἲࢆᥦࡍࡿ㸬࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࡣ㸪ࢹ࣮ࢱࡢ
㸱㸬ᡭ᭩ࡁᩘᏐุูၥ㢟
ศᕸࢆ௬ᐃࡋ࡞࠸ࣀࣥࣃ࣓ࣛࢺࣜࢵࢡ࡞ᡭἲ࡛࠶
ᡭ᭩ࡁᩘᏐุูၥ㢟ᑐࡋ࡚㸪࣮࢝ࢿࣝᐦᗘ᥎
ࡾ㸪」㞧࡞⌧㇟ࢆ☜⋡ⓗࡘᰂ㌾ࡽ࠼࡚⾲⌧
ᐃἲࢆ⏝࠸ࡓ㠀⥺ᙧุูჾ㸦ᩘᏐุูჾ㸧ࢆᵓ⠏
࡛ࡁࡿࡇࡽᵝࠎ࡞ᛂ⏝ࡀᮇᚅࡉࢀࡿ㸬ᮏ✏࡛
ࡋ㸪ࡑࡢุูᛶ⬟ࢆ SVM ẚ㍑ࡋ㸪ホ౯ࡍࡿ㸬ᐇ
ࡣ㸪ᩘᏐุูၥ㢟ᑐࡍࡿ 2 ࡘࡢ≧ἣୗ࡛ࡢ SVM
ࡣ㸪࣮࢜ࣉࣥࢯ࣮ࢫࡢ⤫ィゎᯒࢩࢫࢸ࣒ R[3]
ࡢẚ㍑ࢆ㏻ࡋ࡚ KDE ࡢ≉ᛶࢆ᳨ドࡍࡿ㸬
ࢆ⏝࠸ࡓ㸬
ᐇ㦂⏝ࢹ࣮ࢱ
U.S. Postal Service㸦USPS㸧ZIP code datasets
㸰㸬࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࢆ⏝࠸ࡓ㠀⥺ᙧุูჾ
ከኚ㔞࣮࢝ࢿࣝᐦᗘ᥎ᐃἲ(1)ࡣ㸪ࢹ࣮ࢱⅬ࢝
ࡢᡭ᭩ࡁᩘᏐࢆᑐ㇟ࡋࡓ㸬USPS ࢹ࣮ࢱࡣ 1 ࡘ
࣮ࢿࣝ㛵ᩘࡤࢀࡿᇶᗏ㛵ᩘࢆ㓄⨨ࡋ㸪ࡑࢀࡽ
ࡢᩘᏐࢆ 1616 ࣆࢡࢭࣝࡢࢢ࣮ࣞࢫࢣ್࣮ࣝ㸦-1
ࢆ㡿ᇦෆ࡛㊊ࡋྜࢃࡏࡿࡇࡼࡗ࡚ࡽ࡞᥎
ࡽ 1 ࡢ⠊ᅖࡢ್㸧ࡼࡗ࡚⾲ࡍṇつࡉࢀࡓࢹ
ᐃ㔞ࢆᚓࡿࡶࡢ࡛࠶ࡿ㸬
n
1�
fˆ(x; H) =
KH (x − Xi )
n
࣮ࢱ࡛࠶ࡿ㸬ࢹ࣮ࢱᩘࢆ⾲ 1 グ㍕ࡍࡿ㸬
(1)
i=1
⾲㸯 USPS ࢹ࣮ࢱࢭࢵࢺ
ࢹ࣮ࢱ
ࡓࡔࡋ㸪
=
1, ·· ·, n)
ࡣࢹ࣮ࢱ࣋ࢡࢺࣝ㸪
Xi (i
KH ࡣࣂࣥࢻᖜ⾜ิ H ࢆࡶࡘ࣮࢝ࢿࣝ㛵ᩘ࡛࠶ࡿ㸬
Ꮫ⩦⏝
ከኚ㔞ࢹ࣮ࢱᑐࡍࡿ KDE ࡣ㸪ᵓ⠏ࡢᅔ㞴ࡉࡸ
ࢹ࣮ࢱ
᥎ᐃ⢭ᗘࡢၥ㢟ࡽ┤᥋ⓗ࡞㐺⏝ࡣ㞴ࡋ࠸㸬ࡑࡇ
Ꮫ⩦⏝
࡛㸪㸯ኚ㔞 KDE ࡢ✚࡛⾲⌧ࡍࡿࣉࣟࢲࢡࢺ(✚ᆺ)
ホ౯⏝
ホ౯⏝
0
1194
359
5
556
160
1
1005
264
6
664
170
2
731
198
7
645
147
3
658
166
8
542
166
4
652
200
9
644
177
࣮࢝ࢿࣝ᥎ᐃἲࡼࡗ࡚㏆ఝࡍࡿ㸬ࡓ࠼ࡤ㸪2
ࢢ࣮ࣞࢫࢣ್࣮ࣝࡣ㸪0 ࡽ 1 ࡢ㛫ࡢ್ኚࡋ
ኚ㔞ࣉࣟࢲࢡࢺ KDE ࡣ(2)ᘧ࡛ᐃ⩏ࡉࢀࡿ㸸
☜⋡್ࡋ࡚⏝ࡋࡓ㸬
fˆ(z; hx , hy ) =
1�
Khx (x − Xi )Khy (y − Yj )
n i,j
ࣅ࣮ࣥ࢝ࢿࣝᐦᗘ᥎ᐃἲ
(2)
㏻ᖖࡢ࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࡣ㸪ࢹ࣮ࢱⅬᑐࡋ
࡚ᇶᗏ࡞ࡿ࣮࢝ࢿࣝ㛵ᩘࢆ⏝࠸ࡿ㸬࢝࢘ࣥࢺ(㢖
ᗘ)ࢹ࣮ࢱࡢሙྜࡣ㸪࢝࢘ࣥࢺ್ẚࡍࡿ㔜ࡳࡘ
௨ୗ㸪(2)♧ࡋࡓ KDE ࢆ⏝ࡋࡓ㠀⥺ᙧุ
ࡁࡢ࣮࢝ࢿࣝ㛵ᩘࢆ⏝࠸ࡓࣅ࣮ࣥ࢝ࢿࣝ᥎ᐃἲ
ูࡢὶࢀࢆ♧ࡍ㸬
ࡀ⏝ࡉࢀࡿ㸬USPS ࢹ࣮ࢱࡢᛶ㉁ࡽ㸪1616
1) ࢡࣛࢫูṇゎ☜⋡ศᕸࡢసᡂ
ࣆࢡࢭࣝࢢ࣮ࣞࢫࢣ್࣮ࣝࡢ㧗ࡉࢆᣢࡘࢹ࣮ࢱ
Ꮫ⩦⏝ࢹ࣮ࢱᑐࡋ㸪ࢡࣛࢫࡈ KDE ࢆ㐺⏝ࡋ㸪
ࡳ࡞ࡋ㸪ࣅ࣮ࣥ࢝ࢿࣝᐦᗘ᥎ᐃἲࢆ㐺⏝ࡋࡓ㸬
ࢡࣛࢫู☜⋡ศᕸࡢ᥎ᐃࢆ⾜࠺㸬
ᩘᏐุู࠾ࡅࡿࣅ࣮ࣥ࢝ࢿࣝ᥎ᐃἲࡣ(3)ᘧ
2) ホ౯☜⋡ศᕸࡢసᡂ
࡛ᐃ⩏ࡉࢀࡿ㸬
fˆ(z; hx , hy )
�
=
qij Khx (x − iδ)Khy (y − jδ)
ุูᑐ㇟ࡢホ౯⏝ࢹ࣮ࢱᑐࡋ㸪KDE ࢆ㐺⏝ࡋ㸪
ホ౯☜⋡ศᕸࡢ᥎ᐃࢆ⾜࠺㸬
3) ṇゎศᕸホ౯ศᕸ㛫࡛ࡢ㢮ఝ㔞ࡢ⟬ฟ
㸦3㸧
i,j
�
ࡓࡔࡋ㸪 qij ≡ gij/ i,j
gij ࡛࠶ࡾ㸪 gij ,iδ, jδ
ࡣ㸪
ࡑࢀࡒࢀ (i, j) ࣅࣥ࠾ࡅࡿࢢ࣮ࣞࢫࢣ್࣮ࣝ
ุูᑐ㇟ࢹ࣮ࢱᇶ࡙ࡃホ౯ศᕸྛࢡࣛࢫࡢṇ
ゎศᕸ㛫ࡢ㢮ఝ㔞ࢆᖹᆒ✚ศㄗᕪࡼࡗ࡚⟬
ฟࡍࡿ㸬
ࣅࣥ୰Ⅼ࡛࠶ࡿ㸬
4) ࢡࣛࢫࡢุู
ṇゎ☜⋡ศᕸࡢసᡂ
ホ౯ศᕸ᭱ࡶ㢮ఝࡍࡿṇゎศᕸࡢࢡࣛࢫุู
⾲ 1 ♧ࡋࡓᏛ⩦⏝ࢹ࣮ࢱࢆᩘᏐࡈྛࣅࣥ
ࡍࡿ㸬
࡛ᖹᆒࡋ㸪ࣅ࣮ࣥ࢝ࢿࣝᐦᗘ᥎ᐃࢆ⾜࠺㸬᥎ᐃ
Ⅼࡣྛࣅࣥࡢ୰Ⅼࢆຍ࠼ࡓ 3333 Ⅼ㸪࣮࢝ࢿࣝ
ࣀࢬࢆྵࡴࢹ࣮ࢱࡢุูᛶ⬟
㛵ᩘࡣ㸪2 ኚ㔞ṇつࣉࣟࢲࢡࢺ࣮࢝ࢿࣝ㸪ࣂࣥࢻ
ᖜࡣ㸪ྛḟඖࡶ 1 タᐃࡋࡓ㸬ࢢ࣮ࣞࢫࢣ࣮ࣝ
್ࡀᶆᮏᩘ࡛ࡣ࡞࠸ࡓࡵ㸪⌮ㄽⓗ࡞㆟ㄽࡣ┬ࡃ㸬
.(51(/
᥎ᐃࡋࡓྛᩘᏐࡢṇゎ☜⋡ศᕸࢆᅗ㸯♧ࡍ㸬
690
ࣀࢬࢆྵࡴࣆࢡࢭࣝࡢྜ
ᅗ 2 ࣀࢬࢆྵࡴࢹ࣮ࢱࡢุูᛶ⬟
ᛶ⬟ࡢపୗࡣᑡ࡞࠸ࡇࡀࢃࡿ㸬
ᅗ㸯 ᥎ᐃࡋࡓṇゎ☜⋡ศᕸ
ࢡࣛࢫࢱࣜࣥࢢࡼࡿุูᛶ⬟ࡢᨵၿ
ホ౯⏝ࢹ࣮ࢱࢆ⏝࠸ࡓᩘᏐุู
Ꮫ⩦⏝ࢹ࣮ࢱᛂࡌ࡚ᩘᏐẖ」ᩘಶࡢṇゎ☜⋡
3.3 ࡛సᡂࡋࡓ”0”ࡽ”9”ࡲ࡛ࡢṇゎ☜⋡ศᕸ
ศᕸࢆసᡂࡋ㸪ุูᛶ⬟ࡢᨵၿࢆヨࡳࡓ㸬Ꮫ⩦⏝ࢹ
ホ౯⏝ࢹ࣮ࢱ࡛సᡂࡋࡓホ౯☜⋡ศᕸࡢ㢮ఝ
࣮ࢱࢆᩘᏐࡈ k-means ἲࡼࡗ࡚ࢡࣛࢫࢱࣜ
ᗘࢆ(4)ᘧ࡛࠼ࡽࢀࡿᖹᆒ✚ศㄗᕪ(MISE)
ࣥࢢࡋ㸪ࢡࣛࢫࢱ࣮ࡈṇゎᐦᗘศᕸࢆ᥎ᐃࡍࡿ㸬
ࢆ⏝࠸࡚⟬ฟࡍࡿ㸬
ホ౯⏝ࢹ࣮ࢱ࡛ࡢุูࡢ㝿ࡣ㸪ྛࢡࣛࢫࢱ࣮ࡀᡤ
MISE[fˆ(·; hx , hy )]
� ��
�2 �
ˆ
= E f (z; hx , hy ) − f (z)
dz
ᒓࡋ࡚࠸ࡓᩘᏐุูࡍࡿ㸬ࢡࣛࢫࢱ࣮ᩘࢆኚࡉ
(4)
ࡏࡓࡁࡢุูᛶ⬟ࢆᅗ 3 ♧ࡍ㸬
ࢡࣛࢫࢱ࣮ᩘุูᛶ⬟
MISE ࡢ᭱ࡶᑠࡉࡗࡓ(᭱ࡶ㢮ఝࡋࡓ)ࢡࣛࢫࢆ
ホ౯ࢹ࣮ࢱࡢࢡࣛࢫࡋ࡚ศ㢮ࡍࡿ㸬
690
ุูᛶ⬟ࡢ᳨ド
ᡭ᭩ࡁᩘᏐุูၥ㢟ᑐࡍࡿྛᩘᏐࡢุูᛶ⬟ࢆ
⾲ 2 ♧ࡍ㸬
.(51(/
㸲㸬⤖ᯝ
ุู⋡
࣮࢝ࢿࣝ
690
ࢡࣛࢫࢱ࣮ᩘ
⾲ ホ౯⏝ࢹ࣮ࢱࡢุู⋡
⾲ 2 ♧ࡍࡼ࠺㸪
ᅗ 3 ࢡࣛࢫࢱ࣮ᩘุูᛶ⬟
ᅗ 3 ࡼࡾ㸪ࢡࣛࢫࢱ࣮ᩘࢆቑࡸࡍ㸪SVM ㏆
࣮࢝ࢿࣝ᥎ᐃἲࢆ⏝࠸ࡓ
ᩘᏐุูჾࡣ㸪ุูࡍࡿホ
౯⏝ࢹ࣮ࢱࡢᩘᏐࡼࡗ
࡚㸪ุู⋡ࡤࡽࡘࡁࡀ࠶
㸳㸬⪃ᐹ
ࡿࡶࡢࡢ㸪ᖹᆒࡋ࡚ 0.82
࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࢆ⏝࠸ࡓࣃࢱ࣮ࣥㄆ㆑ᡭἲ
ࡘ࠸࡚㸪ᡭ᭩ࡁᩘᏐุูࡢࢆྲྀࡾୖࡆ᳨࡚ウࡋࡓ㸬
࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࡼࡿᩘᏐุูჾࡣ㸪SVM ࡢ
ุูᛶ⬟ࡼࡾࡶຎࡗ࡚࠸ࡓࡀ㸪ࣀࢬࡢ࠶ࡿࢹ࣮ࢱ
0.94 ࠸࠺㧗ᛶ⬟࡞⤖ᯝ
ᑐࡋ࡚ࡣ㸪SVM ࡼࡾࡶᛶ⬟పୗࡢᙳ㡪ࡣཷࡅ
࡛࠶ࡗࡓ㸬
ィ
࡛࠶ࡗࡓ㸬୍᪉㸪SVM ࡣ
ࡍ࡚ࡢᩘᏐ࠾࠸࡚㸪ุ
ู⋡ 0.9 ௨ୖ㸪ᖹᆒࡋ࡚
㡹ᛶࡢ᳨ド
࠸ุูᛶ⬟ࡀᚓࡽࢀ࡚࠸ࡿࡇࡀࢃࡿ㸬
ࡃࡃ㸪㡹ᛶࡀࡳࡽࢀࡓ㸬ุูࡍࡿᩘᏐࡼࡗ࡚ᛶ
⬟ࡤࡽࡘࡁࡀぢࡽࢀࡿࡓࡵ㸪ࢡࣛࢫࢱࣜࣥࢢࡼ
ࡗ࡚」ᩘಶࡢṇゎᐦᗘศᕸࢆ‽ഛࡍࡿࡇ࡛㸪SVM
࣮࢝ࢿࣝᐦᗘ᥎ᐃἲࡼࡿ㠀⥺ᙧุูჾࡢ㡹
ᛶࢆㄪࡿࡓࡵ㸪ࣀࢬࢆຍࡋࡓᡭ᭩ࡁᩘᏐุู
㛵ࡍࡿᐇ㦂ࢆ⾜ࡗࡓ㸬
ᡭ᭩ࡁᩘᏐࢹ࣮ࢱࣀࢬࢆຍ࠼㸪ุูᛶ⬟ࢆホ
౯ࡍࡿ㸬1616 ࣆࢡࢭࣝࡢ࠺ࡕ㸪ࣀࢬࡢࡁࡉࢆ
ṇつᩘ(ᖹᆒ 0㸪ᶆ‽೫ᕪ 0.5)࡛࠼㸪ࣀࢬࢆ
ຍࡍࡿࣆࢡࢭࣝࡢྜࢆቑࡸࡋࡓࡁࡢุูᛶ⬟ࡢ
ࢢࣛࣇࢆᅗ 2 ♧ࡍ㸬
ᅗ 2 ࡼࡾ㸪SVM ࡢุูᛶ⬟ࡣࣀࢬࡢ㔞ࡶ
పୗࡍࡿࡀ㸪࣮࢝ࢿࣝ᥎ᐃἲࡼࡿᩘᏐุูჾࡢ
ࢃࡽ࡞࠸ุูᛶ⬟ࡀᚓࡽࢀࡓ㸬ᐇ⏝ྥࡅ㸪
ࡉࡽ࡞ࡿᨵၿࡀᮇᚅ࡛ࡁࡿ㸬
ཧ⪃ᩥ⊩
[1] V.Vapnik, "The Nature of Statistical Learning Theory",
Springer, (1995).
[2] M. P. Wand, M. C. Jones "Kernel Smoothing", Monographs
on Statistics and Applied Probability, Chapman & Hall, (1995)
[3]R Development Core Team “R: A language and
environment for statistical computing”, R Foundation for
Statistical
Computing,
Vienna,
Austria.
URL
http://www.R-project.org/.
Asymptotics for M-Estimators in Time Series
Yan LIU
Department of Pure and Applied Mathematics, Waseda University, Tokyo 169-8555
November 28, 2013
at Kanazawa University
In this talk, we give conditions which guarantee the asymptotic normality of M-estimator
based on the observations from autoregressive (AR) models in the time series analysis by
minimizing some convex objective functions. The objective function is defined by twovariate function since the scale parameter of time series model has to be estimated as well
as the coefficient parameters. In addition, we do not assume the differentiability around
the true parameter of the objective functions. The results are extended to Mm estimates.
Examples of the results are also provided.
First, we reviewed the history on the development of M-estimation from the point of view
of relaxing the condition of the differentiability around the true value. It has been a half
of a century since Hodges and Lehmann (1963) proposed a Lemma on the asymptotic normality for the estimation of location before Huber, who changed the idea into the concept
of M-estimation. More precisely, Huber showed the asymptotic normality of M-estimators
for convex real-valued function and investigated the robustness of the estimators in 1964.
Koenker and Bassett (1978) provided the theory on the regression quantiles by employing
the idea of Huber. Pollard (1991) reinvestigated the properties of convex objective function and proposed the convexity lemma. Niemiro (1992), following the line in Habermann
(1989), established a Bahadur-type representation of the estimators. Further, Bose (1998)
generalized the results in Niemiro (1992) to the case of Mm estimators including Oja’s median, univariate robust scale estimator of Bickel and Lehmann and Hodges and Lehmann’s
estimator.
To evaluate the efficiency and robustness of the estimators, we focused on Martin and
Yohai (1985). They summarized the concepts into four types such as Tukey’s efficiency
robustness, Huber’s min-max robustness, Hampel’s qualitative robustness and Tukey’s resistance. We did not enter into the last two concepts here. Tukey’s efficiency robustness is
defined as the ratio of the Cramer-Rao lower bound at the parametrized joint distribution
of the samples to the variance of the statistic derived from the samples. In the asymptotic
case, the ratio can be written as one over the product of the asymptotic Fisher information and the asymptotic variance of the statistic. Huber (1981) formalized the min-max
robustness on the asymptotic variance of the statistic at the distributions of the samples
and call the statistic min-max robust if the statistic is the infimum over all statistics of
the supremum over all distributions of the samples of the asymptotic variance. Since the
1
underlying distribution is not known generally, Huber (1964) and Huber (1973) proposed
minimizing the objective function of i.i.d. samples normalized by its scale parameter.
However, the scale parameter is also unknown in time series analysis. We further generalized M-estimators to be a two-variate function, estimating the scale parameter and
coefficient parameters of AR models simultaneously. The usage of two-variate function corresponds to the Whittle’s estimators, which is optimal in higher order efficiency in the
frequency domain in time series analysis. We suppose the objective function to be convex
but may not be differentiable around the true parameter. In this case, the theory contains
the L1 norm, which means that the quantile estimators such as least absolute deviation
(LAD) estimator are also considered. Thus this class of M-estimators is larger than the
Whittle’s estimator since the theory in the frequency domain mainly depends on L2 norm.
The estimators may be more robust beyond those estimators whose objective function
defined in L2 norm if the underlying distribution of the samples has a sufficiently large
deviation. The approach for L1 norm, which receives broad attention, has been developed
in the i.i.d. case, regression case and also recently in time series case (see Koenker (2005)).
The conditions in this talk for asymptotic normality of M-estimators in time series are
chiefly derived from both Hodges and Lehmann (1963)’s and Niemiro (1992)’s conditions.
The former ones are given by Hodges and Lehmann and rearranged in Inagaki and Kondo
(1980). The conditions are not directly assumed on the objective function but on the score
function. The convergence in probability of the difference between the score function with
√
true parameter and its n-neighborhood is difficult to show, and we employ the latter one
for proof. Although Niemiro (1992) is oriented to show that M-estimators have Bahadurtype representation, they gave a very nice structure of proof for asymptotic normality of
the estimators with convex objective functions. When the objective function is convex,
there exists a subgradient function of the objective funtion. The subgradient function can
be regarded as the score function. The asymptotic normality of M-estimators is derived
from the properties of the subgradient function. The main difficulty in time series case
is that only L2 integrability of the score function is not enough. To show the asymptotic
normality, we have further to guarantee that the covariance between the components of
the score function is also summable, like in Billingsley (1968). Under the conditions, we
obtained the desirable results without the differentiability around the true parameter. We
also derived the asymptotic variance of the estimators assuming the differentiability of the
estimators. Our results agree with the existing results, which were given in the talk as
examples. The extension of the results to Mm estimators was also mentioned in the talk,
which can be shown by the properties of U-statistics given in Yoshihara (1976) or Denker
and Keller (1983). The idea was generalized from Hoeffding (1948) and Hoeffding (1961).
The contribution of the work is that we generalized the approach of M-estimation to
estimate the scale parameter and the coefficient parameters simultaneously and defined the
objective function convex, so that the objective function do not have to be differentiable
around the true parameters. We also gave the conditions for asymptotic normality of the
estimators, which in turn can be generalized to the approach in the frequency domain.
2
Statistical Inference Associated with the Fractional Brownian Motion
and Related Processes
田中 勝人（学習院大学経済学部）
1. フラクショナル O-U 過程
次の連続時間確率過程を考察の対象とした．
dYH (t) = α YH (t) dt + dBH (t),
(α ≤ 0),
(t ∈ [0, M ])
ここで，{BH (t)} は，Hurst パラメータ H (1/2 ≤ H < 1) の fBm (フラクショナル Brown 運動)
であり，{YH (t)} は，fO-U (フラクショナル Ornstein-Uhlenbeck) 過程と呼ばれる．ここでは，ドリ
フト・パラメータ α (≤ 0) の推定問題と検定問題を考察した．
2. 推定
推定量としては，LSE (最小 2 乗推定量)，MLE (最尤推定量)，MCE (最小コントラスト推定量)
を取り上げて，それぞれの推定量の分布を考察した．LSE については，次の 3 つの推定量を取り上
げた．
!M
!M
YH (t) dBH (t)
0 YH (t) dYH (t)
= α+ 0 !M
α̂1 =
!M 2
YH (t) dt
YH2 (t) dt
0
0
!M
!M
YH (t) δYH (t)
YH (t) δBH (t)
0
0
=
α
+
α̂2 =
!M 2
!M 2
YH (t) dt
YH (t) dt
0
0
"
$−1/2H
# M
1
α̂3 = −
YH2 (t) dt
HΓ(2H)M 0
ここで，α̂1 の定義の分子で使われている積分は，Wick 積分であり，α̂2 の積分では，通常のリーマ
ン和による定義が使われている．この中で，α̂1 と α̂3 は，観測時間 M が無限大のとき一致性をもつ
が，α̂2 は一致性をもたない．また，α̂3 の方が，α̂1 よりも効率的であることが知られている．なお，
有限の M に対して，これらの分布関数や密度関数の計算は，今のところ，未解決な問題である．
MLE α̃MLE は，フラクショナルな場合の Girsanov の定理により，
!M
QH (t) dZH (t)
UH
=
α̃MLE = !0M
,
2
VH
Q (t) dvH (t)
0
H
で定義される．ここで，3 つの確率過程 ZH (t), QH (t), vH (t) が新たに導入されている． MLE の分
布を計算するために，UH と VH の同時積率母関数 m(θ1 , θ2 ) を，Klepstyna and Le Breton (2002)
の結果を使って導出して，分布関数は次の公式に従って数値計算することができた．
#
1 ∞ 1
1
Im [m(−iθ, iθx)] dθ
P (α̃MLE < x) = +
2 π 0 θ
ただし，
m(θ1 , θ2 )
=
%&
'
α + θ1
µM
(α + θ1 )2
−
sinh µM
1+
cosh2
µ2
2
µ
&
&
(
'
'
µM
µM
πM
(α + θ1 )2
+
I−H
−
IH−1
4 sin πH
µ
2
2
&
&
'
')*−1/2
µM
µM
+µI1−H
,
IH
2
2
e−M(α+θ1 )/2
√
ここで，µ = α2 − 2θ2 であり，Iν (z) は，階数 ν の変形 Bessel 関数である．
他方，MCE α̂MCE は，
−M/2
α̂MCE = ! M
Q2H (t) dvH (t)
0
で定義される．MLE と MCE については，密度関数を数値積分により計算して，グラフ表示した．
また，M → ∞ のときの漸近分布については，次の結果を示した．
√
√
M (α̃ − α) & M (α̃MCE − α) → N(0, −2α)
(α < 0)
なお，α = 0 の場合は，次の結果が成り立つ．
D
M α̂MCE (H, M ) = α̂MCE (H, 1)
(α = 0)
この結果は，他の推定量についても成り立つ．したがって，α = 0 の場合，推定量は，M − consistent
であるが，漸近分布は非正規である．
3. 検定
パラメータ α に関する検定問題
H0 : α = 0
vs.
H1 : α < 0
を考察した．この問題は，離散的な場合の拡張された単位根近接モデル
yj = ρyj−1 + vj ,
ρ=1+
α
,
T
(1 − L)H−1/2 vj = εj ,
(j = 1, · · · , T ),
における単位根検定の連続時間モデル・バージョンであると解釈することができる．
検定は，MLE と MCE に基づく統計量を考えて，検出力を計算した．そのために，まず，帰無
分布の分位点を計算して，次の結果を得た．
Probability
H
0.01
0.05
0.1
of
a
smaller
value
0.5
0.9
0.95
0.99
−0.853
−0.836
−0.636
0.928
0.899
0.767
1.285
1.250
1.084
2.033
1.989
1.791
−1.721
−1.606
−1.027
−0.418
−0.378
−0.209
−0.302
−0.272
−0.149
−0.179
−0.161
−0.087
α̃MLE
0.5
0.7
0.9
−13.696
−13.608
−12.988
−8.039
−7.964
−7.415
−5.714
−5.648
−5.154
α̂MCE
0.5
0.7
0.9
−14.510
−14.314
−13.196
−8.856
−8.673
−7.630
−6.533
−6.359
−5.376
そして，検出力を計算することにより， MLE に基づく検定の方が，検出力が高くなることを見
出した．
4. 今後の課題
ここで考察した fO-U 過程は，通常の O-U 過程を拡張したものであり，長期記憶性をもつよう
な確率過程を描写するモデルとして，今後，ますます研究が進められると思われる．現時点で未解決
な問題は，パラメータ α の最小 2 乗推定量の分布の導出，離散時間モデルの推定との関係，Hurst
パラメータ H の推定問題などである．これらについても，今後の研究で明らかにしていきたい．
1
1
2
2
1.
t
β(t)
Satoh and Yanagihara (2010, Amer. J. Math. Management Sci. 30:243-256)
t∈R
(2013,
t ∈ [a, b]
42:11-21)
β(t) = θ ! x(t)
2.
θ̂ ∼ Nq (θ, Ω)
θ̂
!
β(t)
!
!
λ(t) = x(t) Ωx(t)
2
β̂(t) = θ̂ x(t) ∼ N(β(t), λ (t))
Ω̂
β(t)
θ = (θ1 , . . . , θq )!
100(1 − α)%
"
#
I1−α (t|uα ) = β̂(t) − uα λ̂(t), β̂(t) + uα λ̂(t) ,
I1−α (t|uα )
uα
1−α
β(t)
Pr(T ≤ uα ) = 1 − α,
[a, b]
Pr (β(t) ∈ I1−α (t|uα )) =
T = sup |T (t)|,
t∈[a,b]
T (t) =
β̂(t) − β(t)
,
λ(t)
uα
K
T (t)
$
x(t)! Ω̂x(t)
Pr (η(t) ∈ I1−α (t|uα ))
t ∈ [a, b]
uα
λ̂(t) =
Ω
Tk = T (tk ) (k = 1, . . . , K)
T ≈ Tmax =
K
Tmax
IK (a, b) = {tk , k = 1, · · · , K|a = t1 < t2 < · · · < tK = b}
T
max |T (t)| = max {|T1 |, |T2 |, . . . , |TK |} .
t∈IK (a,b)
T
Pr(Tmax ≤ uα ) = 1 − α ⇔ Pr (|T1 | ≤ uα , . . . , |TK | ≤ uα ) = 1 − α
t = (T1 , . . . , TK )!
t = (T1 , . . . , TK )! ! NK (0, R)
R = DCΩC ! D, C = (x(t1 ), . . . , x(tK ))! , D = diag(λ(t1 ), . . . , λ(tK ))−1
uα
ΦK (x|R)
K
ΦK (x|R) =
ΦK (x|R) = Pr(−x ≤ t1 ≤ x, . . . , −x ≤ tK ≤ x|R)
1
1/2
|R| (2π)K/2
%
x
−x
1
···
%
x
−x
1 !
e− 2 t R
−1
t
dt,
1
ΦK (xα |R) = 1 − α
xα
Pr (β(t) ∈ I1−α (t|xα )) ≈ 1 − α
Genz and Bretz (2009)
Genz and Bretz (2002)
3.
I1−α (t|xα )
R
Bretz and Genz 2001)
t ∈ [a, b]
uα = xα
Genz (1992, 1993)
mvtnorm
(Hothorn,
1000
Satoh and Yanagihara (2010)
(2) t ∈ R
(1)
(
(Satoh and Yanagihara 2010) (3) t ∈ [a, b]
2013)
and Roy (1964)
Potthoff
16
2
4
)
.
11
(mm)
t
y(t)
a = 1)
E[y(t)]
!
E[y(t)] = a β(t) = β1 (t) + β2 (t)a,
βj (t) =
(8
a
14
(
a=0
β(t) = (β1 (t), β2 (t))
θj! x(t),
!
a = (1, a)!
β1 (t)
β2 (t)
2 !
x(t) = (1, t, t )
(
Laird and Ware 1982)
8
t ∈ [8, 14]
14
0.95
Satoh and Yanagihara (2010)
t
I0.95 (t|2.365)
I0.95 (t|1.960)
t∈R
x0.05 = 2.365
√
c3,0.05 = 2.795
q0.05/2 = 1.960
Satoh and Yanagihara (2010)
1
1
β1 (t)
I0.95 (t|2.795)
Satoh and Yanagihara (2010)
2
2
⤒ࢸ࢟ࢫࢺࢹ࣮ࢱᑐࡍࡿከḟඖᑻᗘἲࡢᛂ⏝
ศᏛ ἨᚿὠᜨࠊᗈᓥᏛ బ⸨୍
㸯ࠋࡣࡌࡵ
㏆ᖺࠊὀ┠ࡉࢀ࡚࠸ࡿࣅࢵࢢࢹ࣮ࢱࡢ୍ࡘ⤒ⓗほ ࡉࢀࡓࢸ࢟ࢫࢺࢹ࣮ࢱࡀ࠶ࡿࠋࢸ
࢟ࢫࢺࢹ࣮ࢱྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡢ⣔ิⓗ࡞ฟ⌧ࣃࢱ࣮ࣥࢆᢳฟࡋࠊࡑࡢศ㢮ࢆ⾜࠺ࡇ
࡛ࠊࢸ࢟ࢫࢺࢹ࣮ࢱࡢ≉ᚩⓗ࡞ኚࡋ࡚せ⣙࡛ࡁࡿྍ⬟ᛶࡀ࠶ࡿࠋ
୍᪉ࠊᶓ㍈ ᐃ㛫ࠊ⦪㍈࣮࣮࢟࣡ࢻฟ⌧ࡢ᭷↓ࢆࡗࡓᩓᕸᅗࢆ⪃࠼ࡿࠊฟ⌧ࣃࢱ
࣮ࣥࡢᢳฟ᪉ἲࡋ್࡚㞳ᩓࢹ࣮ࢱᑐࡍࡿᖹࡀ⪃࠼ࡽࢀࡿࠋᖹࡢᡭἲࡋ࡚ࡣ
ࡁࡃศࡅ࡚ࠊ㸯㸧ᒁᡤ⥺ᙧᅇᖐࢆ⧞ࡾ㏉ࡋ⾜࠺࣮࢝ࢿࣝᖹࠊ㸰㸧ࢫࣉࣛࣥᇶᗏ࡞ࢆ⏝
࠸࡚ᇦⓗ᥎ᐃࢆ⾜࠺ࣃ࣓ࣛࢺࣜࢵࢡᅇᖐࠊࡀ࠶ࡿࠋࡇࡇ࡛ࡣࠊᚋ⪅ูࡉࢀࡿΰྜຠᯝ
ࣔࢹࣝࢆ⏝࠸ࡓ᥎ᐃ᪉ἲ㸦࠼ࡤࠊ Brumback et al. (JASA, 1999)ࠊ బ⸨࣭⏣(ᛂ⏝⤫ィᏛ,
2013)㸧ࢆ㐺⏝ࡋࠊࡑࡢ᥎ᐃ⤖ᯝࢆ⏝ࡋࡓศ㢮᪉ἲどぬࢆᥦࡍࡿࠋ
㸰㸬⤒ⓗࢸ࢟ࢫࢺࢹ࣮ࢱࡢฟ⌧ࣃࢱ࣮ࣥࡢᢳฟศ㢮
ࡲࡎࠊn Ⅼ࠾࠸࡚㞟ࡵࡽࢀࡓࢸ࢟ࢫࢺࢹ࣮ࢱࡽᙧែ⣲ゎᯒࡼࡾ m ಶࡢ㧗ฟ⌧㢖ᗘࡢ
༢ㄒࢆ࣮࣮࢟࣡ࢻࡋ࡚ᢳฟࡍࡿࠋࡑࡋ࡚ࠊᢳฟࡋࡓ m ಶࡢ࣮࣮࢟࣡ࢻࡢྛⅬ࠾ࡅࡿฟ⌧
ࡢ᭷↓ࡽ n ⾜ m ิࡢ 2 ್⾜ิࢆ⏕ᡂࡍࡿࠋḟࠊฟ⌧ࡢ᭷↓ࢆ┠ⓗኚᩘࡋ࡚ࠊ ᐃⅬ
ࢆ⏝࠸ࡓ┤⥺࠸ࡃࡘࡢ⠇Ⅼࢆᣢࡘᢡࢀ⥺ࢆᇶᗏࡍࡿࢭ࣑ࣃ࣓ࣛࢺࣜࢵࢡ࡞ࣟࢪࢫࢸ
ࢵࢡᅇᖐࢆ⾜࠺ࠋࡇࡢࡁࠊᢡࢀ⥺ࣜࢵࢪᆺࡢ⨩๎ࢆ࠼ࡿࡇࡣࠊBrumback et al. (1999)
ࡢᥦᡭἲࡽࠊᢡࢀ⥺ࡢᅇᖐಀᩘࣛࣥࢲ࣒ຠᯝࢆ௬ᐃࡋࡓ୍⯡⥺ᙧΰྜࣔࢹࣝࢆ⪃࠼ࡿ
ࡇ➼ࡋࡃࠊ⨩๎ࣃ࣓࣮ࣛࢱᅇᖐಀᩘࡣྠ᥎ᐃࡉࢀࡿࠋ
ḟࠊ࣮࣮࢟࣡ࢻࡈ᥎ᐃࡉࢀࡓ⌮ㄽ᭤⥺ࡢศ㢮ࢆ⪃࠼ࡿࠋ⌮ㄽ᭤⥺ࢆせ⣙ࡍࡿࡶࡗࡶ
⮬↛࡞せ⣙㔞ࡣᅇᖐಀᩘ࣋ࢡࢺ࡛ࣝ࠶ࡿࠋࡋࡋ࡞ࡀࡽࠊᅇᖐಀᩘ࣋ࢡࢺࣝࡢྲྀࡾᚓࡿ⠊ᅖࡣ
ᗈࡃࠊࡲࡓ⌮ㄽ᭤⥺ࡀ㢮ఝࡋ࡚࠸࡚ࡶᅇᖐಀᩘ࣋ࢡࢺࣝࡋ࡚ࡣࡁ࡞ᕪࢆᣢࡘࡇࡀ࠶ࡿࠋ
ࡑࡇ࡛ࠊᮏ✏࡛ࡣ⌮ㄽ᭤⥺ࡢ್ࡀ(0, 1)༊㛫㝈ᐃࡉࢀࡿࡇ╔┠ࡋࠊ⠇Ⅼࡢ⨨ ᐃⅬ
ࡢ୧➃࠾ࡅࡿ⌮ㄽ್ࢆ⌮ㄽ᭤⥺ࡢせ⣙㔞ࡋ࡚⏝࠸ࡿࡇࢆ⪃࠼ࡿࠋࡇࡢࡁࠊせ⣙㔞ࡢḟ
ඖࡣᅇᖐಀᩘ࣋ࢡࢺࣝࡢ㛗ࡉ➼ࡋࡃࠊฟ⌧ࣃࢱ࣮ࣥࡀ㢮ఝࡋ࡚࠸ࢀࡤせ⣙㔞ࡶ㢮ఝࡍࡿࠋ
ࡉࡽࠊせ⣙㔞ᑐࡋ࡚ k-means ἲࢆ㐺⏝ࡍࡿࡇ࡛ฟ⌧ࣃࢱ࣮ࣥࢆศ㢮ࡍࡿࠋ⤖ᯝࡋ࡚
ᚓࡽࢀࡿ௦⾲ⓗ࡞せ⣙㔞ࡣ⤒ⓗࣉࣟࢵࢺࡍࡿࡇ࡛௦⾲ⓗ࡞ฟ⌧ࣃࢱ࣮ࣥࡋ࡚どぬ
࡛ࡁࡿࠋࡲࡓࠊฟ⌧ࣃࢱ࣮ࣥࡢᙧ≧ࡣኻࢃࢀࡿࡀࠊࡍ࡚ࡢ࣮࣮࢟࣡ࢻࡢせ⣙㔞ࡑࡢ௦⾲Ⅼ
ᑐࡋ࡚ከḟඖᑻᗘἲࢆ⾜࠺ࡇ࡛ศ㢮⤖ᯝࡣ➃ⓗ࡞ 2 ḟඖᩓᕸᅗࡋ࡚⾲⌧࡛ࡁࡿࠋ
㸱㸬ᐇࢹ࣮ࢱࡢᛂ⏝
ᥦࡋࡓ⤒ⓗほ ࡉࢀࡓࢸ࢟ࢫࢺࢹ࣮ࢱࡢฟ⌧ഴྥࡢᢳฟศ㢮ࡢ᪉ἲ࠾ࡼࡧどぬ
ࡘ࠸࡚ࠊᐇࢹ࣮ࢱ㐺⏝ࡋࡓࢆ⤂ࡍࡿࠋࡇࡇ࡛ࡣࠊᗈᓥᕷ㛗ࡀୡ⏺ྥࡅ࡚Ⓨ⾲ࡋࡓᖹ
ᐉゝࡢෆࠊ1950 ᖺࢆ㝖ࡃ 1947 ᖺ-2012 ᖺࡢ 65 Ⅼほ ࡉࢀࡓ᪥ᮏㄒࢸ࢟ࢫࢺࢹ࣮ࢱࢆ⏝
࠸ࡿࠋࡇࡢࢸ࢟ࢫࢺࢹ࣮ࢱ୰ࡢ࣮࣮࢟࣡ࢻࡢฟ⌧ࡢ᭷↓ࡀ㛫ࡶࡢࡼ࠺ኚࡍࡿࡢ
ࠊࡑࡋ࡚ࠊࡢ࣮࣮࢟࣡ࢻࡀඹ㏻ࡢฟ⌧ഴྥࢆᣢࡘࡢ㛵ᚰࡀ࠶ࡿࠋ⤒ࢸ࢟ࢫࢺࢹ࣮ࢱ
RMeCab㸦▼⏣, 2008㸧ࢆ⏝࠸࡚ᙧែ⣲ゎᯒࢆ⾜࠸ࠊ㧗ฟ⌧㢖ᗘࡢୖ ༢ㄒࢆᢳฟࡋࡓࠋ
ࡇࡢ⤖ᯝࠊᚓࡽࢀࡓほ ್ࡽ࡞ࡿ ⾜ ิࡢ ್⾜ิᑐࡋ࡚ࠊ⤫ィゎᯒࢯࣇࢺ࢙࢘ 5
㸦5&RUH7HDP㸧ࢆ⏝࠸࡚ࠊᥦ᪉ἲࡼࡿ࣮࣮࢟࣡ࢻ㛵ࡍࡿせ⣙ࢆ⾜ࡗࡓࠋ࣮࢟࣡
࣮ࢻࡈࡢฟ⌧ഴྥࢆ ࡘࡢࢡࣛࢫࢱศ㢮ࡍࡿࠊࢡࣛࢫࢱ࣮ྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡑࡢ
ᩘࡣࠊฟ⌧☜⋡ࡀቑຍ࣭ῶᑡࡍࡿ &㹹ᙉ࠸ࠊ㌷⦰ࠊᐇ㦂ࠊᅜᐙࠊ㐃ᖏࠊ⏕Ꮡ㹻ࡢ ಶࠊ ๓
ᚋ࡛ࡰ୍ᐃࡢ☜⋡ࢆಖࡘ &㹹ୡ⏺ࠊᖹࠊே㢮ࠊᗈᓥࠊᡓதࠊᕷẸࠊ㏄࠼ࡿ㹻ࡢ ಶࠊ☜⋡
ࡀቑຍࡋ⥆ࡅࡿ &㹹㒔ᕷࠊᐇ⌧ࠊ㆟ࠊᨻᗓࠊồࡵࡿࠊᕷ㛗ࠊᅜ㐃ࠊᆅ⌫ࠊ᮲⣙ࠊ㛤ദࠊᚰࠊ
ᚚ㟋ࠊㆤࠊኌࠊᮍ᮶ࠊᛮ࠸ࠊୡ⣖㹻ࡢ ಶࠊ ๓ᚋ࡛ࡰ୍ᐃࡢ☜⋡ࢆಖࡘ &㹹ேࠎࠊ
ே㛫ࠊయ㦂ࠊ⚗Ṇࠊດຊࠊ㐨ࠊ᪂ࡓࠊㄋ࠺ࠊỴពࠊグᛕࠊ☜❧ࠊ◚ቯࠊᜏஂ㹻ࡢ ಶࠊ
ᖺ௦ࡢึᮇࡢ㡭ቑຍࡋ୍ᐃࡢ☜⋡ࢆಖࡘ &㹹᰾රჾࠊ⿕⇿ࠊ᰾ࠊᗫ⤯ࠊཎ⇿ࠊࣄࣟࢩ࣐ࠊᅜ
㝿ࠊッ࠼ࡿࠊ≛≅㹻ࡢ ಶ࡞ࡗࡓࠋከḟඖᑻᗘἲࡼࡗ࡚ฟ⌧ഴྥࡢせ⣙㔞ࢡࣛࢫࢱᖹᆒ
ࢆ ḟඖᩓᕸᅗ♧ࡍࡇࡀ࡛ࡁࡿࠋࡇࡢᅗ࡛ࡣࠊᥦ᪉ἲࡼࡿ࣮࣮࢟࣡ࢻࡢฟ⌧☜⋡ࡢ⤒
ⓗ࡞ฟ⌧ഴྥࡢศ㢮⤖ᯝ࣮࣮࢟࣡ࢻࡢ㛵ಀࡀ⾲ࡉࢀ࡚࠾ࡾࠊᚑ᮶ࡢ㛫㡰ᗎ౫Ꮡࡋ࡞
࠸ ್ࡢྠ⏕㉳ᇶ࡙ࡃከḟඖᑻᗘἲࡢ⤖ᯝ␗࡞ࡿⅬࡀከ࠸ࠋ
ࡇࡇ࡛ࡣᖹᐉゝࢆせ⣙ࡍࡿ 5 ࡘࡢࢡࣛࢫࢱ࣮ࢆᐇᏛⓗ࡞どⅬࡽゎ㔘ࡋࠊ≉ᚩ࡙ࡅࢆヨࡳ
ࡿࠋศ㢮 C1 ྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡣࠊᅜ㝿㛵ಀࢆ⾲ࡍ⪃࠼ࡿࠋᮾす෭ᡓࡢ⤊⤖ࡶ࡞࠸ࠊ
ᪧࢯ㐃㸦࢝ࢨࣇࢫࢱࣥ㸧ࢆྵࡴᅜࡀ 1996 ᖺࢆ᭱ᚋࡑࢀࡒࢀᆅୗ᰾ᐇ㦂ࢆṆࡋࡓࠋ
᰾ᐇ㦂ᑐࡍࡿᠱᛕࡸ㌷⦰ᑐࡍࡿᮇᚅࡀ⾲ࢀࡓ࣮࣮࢟࣡ࢻࡢฟ⌧☜⋡ࡣᖺࡶቑຍࡋ
࡚࠸ࡁࠊ1990 ᖺ㡭ࡽῶᑡࡋጞࡵࡿࠋḟࠊศ㢮 C2 ྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡣࠊᗈᓥࡢᬑ㐢ⓗ
࡞ᖹほࢆ⾲ࡍ⪃࠼ࡿࠋほᐹᮇ㛫ࡢ㛫ᕷ㛗ࡀ 8 ᅇ௦ࡋࡓࡶࡢࡢࠊ࣮࣮࢟࣡ࢻࡢฟ⌧☜⋡
ࡣ 0.8 ࡢ㧗࠸್ಖࡓࢀ࡚࠸ࡿࠋḟࠊศ㢮 C3 ྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡣࠊᖹࡢᢸ࠸ᡭࡢ
ኚࢆ⾲ࡍ⪃࠼ࡿࠋ1982 ᖺ㛤ദࡉࢀࡓ➨ 2 ᅇᅜ㐃㌷⦰≉ู⥲ࢆዎᶵᅜࡽ㒔ᕷ࠸࠺
ࡼࡾᑠࡉ࡞ᕷẸ♫ࡢ༢࠾࠸࡚ᖹࡢᐇ⌧ࢆࡵࡊࡍࡇࢆ⾲ࡍ࣮࣮࢟࣡ࢻࡢฟ⌧☜⋡ࡣ
ᖺࠎቑຍࡋ࡚࠸ࡃࠋḟࠊศ㢮 C4 ྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡢฟ⌧☜⋡ࡣほᐹᮇ㛫ࢆ㏻ࡋ࡚ 0.5
⛬ᗘಖࡓࢀ࡚ࡇࡽࠊศ㢮 C2 ࡢࡼ࠺㧗㢖ᗘ࡛ࡣ࡞࠸ࡶࡢࡢࠊࡸࡣࡾࠊ࠶ࡿ⛬ᗘࡢᬑ㐢
ⓗ࡞せ⣲⪃࠼ࡿࠋ᭱ᚋࠊศ㢮 C5 ྵࡲࢀࡿ࣮࣮࢟࣡ࢻࡣࠊ1945 ᖺࡽࡢ✵ⓑࡢ 10 ᖺࢆ
⾲ࡍ⪃࠼ࡿࠋࡇࡢࡼ࠺ᥦ᪉ἲࢆ⏝࠸ࡿࡇ࡛ᖹᐉゝࢆ➃ⓗ 5 ࡘࡢせ⣲ࡋ࡚せ⣙࡛
ࡁࠊ♫⫼ᬒࢆ⪃៖ࡋࡓゎ㔘ࡶྍ⬟࡞ࡗࡓࠋ
㸲㸬ࡲࡵ
⤒ⓗほ ࡉࢀࡓࢸ࢟ࢫࢺࢹ࣮ࢱ࠾࠸࡚ࠊ࣮࣮࢟࣡ࢻࡢฟ⌧ࡢ᭷↓╔┠ࡋࠊ࣮࣮࢟࣡
ࢻࡢ⤒ⓗ࡞ฟ⌧ഴྥࢆΰྜຠᯝࣔࢹࣝࡢࡶ࡛᥎ᐃ᭤⥺ࡋ࡚せ⣙ࡋࠊࡇࢀࢆศ㢮࣭どぬ
ࡍࡿ᪉ἲࢆᥦࡋࡓࠋᐇࢹ࣮ࢱ㐺⏝ࡋ࡚ᚓࡽࢀࡓ⤖ᯝࡢᐇᏛⓗ࡞ゎ㔘ࡣࠊᮏ᪉ἲࡢጇᙜᛶࢆ
᳨ウࡍࡿୖ࡛㔜せ࡞せ⣲࡞ࡾᚓࡓࠋ
㐃⤡ඛ: Ἠ ᚿὠᜨ (Email: [email protected])
E-mail: [email protected]
1. Introduction
In the field of biology theories of aging are roughly divided into two major groups. One is consisting of damage
theories and the other is consisting of program theories. According to damage theories we age because our
systems break down over time. Meanwhile, in the program theories, it is considered that we age because there is
an inbuilt mechanism that tells us to die. If the damage theories are true, we can survive any longer by avoiding
damaging our organism. If the program theories are true, on the other hand, we cannot survive longer than the
upper limit of longevity with any effort.
Considering the above our aim in this study is to see if the damage theories are true, which means that
there exist a upper limit of human longevity, or the programing theories are true, which means that there does not
exist such a limit, by applying the extreme value theory to an analysis of data for Japanese oldest old survivors
obtained using the method of extinct cohort (Wilmoth, Andreev, Jdanov, and Glei, 2007).
2. Generalized Pareto model
Let X 1 , X 2 , , X n be a sequence of independent lifetime of individuals with common distribution function
F and Denote an arbitrary term in the X i sequence by X . Then, according to the extreme value theory, if it
can be considered that u ! 100 is a large enough constant, the distribution function of X u conditional on
X ! u is approximately the generalized Pareto distribution, which is written as
1/ J
J z0
°1 (1 J y / a )
F ( y ; J , a) ®
J 0.
°̄1 exp( y / a )
(2)
An important property of this distribution is that the upper limit of distribution varies depending
on parameter values. That is, if J 0 then the upper limit of the distribution is finite, that is
0 Y a / J { Z , and if J t 0 , then the upper limit of the distribution is infinite, that is, 0 Y f .
Considering those properties, this study is laid out as follows (1) to test the hypothesis J 0 , (2) to estimate the
upper limit of life for each cohort.
Let M ij indicate the numbers of survivors who are 100 i years old ( i
1, 2, , I ) in the year
1976 j ( j 1, 2, , 26 ), where I is supposed to be sufficiently large so that M Ij
numbers of deaths in a year are obtained as difference of M ij , that is, mij
individuals
pj
m
die
in
the
year
1976 j
at
the
age
of
100 i
0 for all j , and the
M ij M i 1, j 1 . Now assume that
with
probabilities
given
by
F (i 1, J j , a j ) F (i, J j , a j ) . Then the contribution of mij to the likelihood is M ij mij qij ij (1 q ij )
M ij mij
,
where qij
p j / (1 F (i, J j , a j )) , hence the likelihood is proportional to
I ,J
i 1, j
§ M ij · mij
M m
¨¨
¸¸ qij (1 q ij ) ij ij .
1 © mij ¹
(2)
Thus models assumed on pij thus fall into the class of binomial regression models, and the parameters in the
generalized Prato distribution are estimated by maximizing the likelihood function (2).
3. Numbers of survivors obtained using the method of extinct cohort
The numbers of deaths by year when they were born are indicated in Vital Statistics in Japan. So, the population
size for a cohort at a certain age is estimated by summing all future deaths for the cohort with the method of
extinct cohort, which is one whose members are assumed to have all died by the end of the observation period.
(U )
The method of extinct cohort is explained as follows. Let Dij
indicate the number of deaths who are
( L)
at the age of i in the year j and were born in the year j i 1 and Dij
indicate the number of deaths
who are at the age of i in the year j and were born in the year j i . Then, assuming that there is no
migration or error, the number of survivors who have reached the age of i , say Sij , is
Sij
Di(,Lj )
f
¦ Di(Uk), j k Di(Lk)1, j k .
k 1
(3)
4. Result
Z obtained by applying the method of extinct
cohort to the data for male and female deaths. As shown in the table, the estimates of J are negative for all male
Table 1 indicates the ML estimates of the parameter J and
cohorts with sufficiently small p-values except the cohort 1968 for male and the cohorts 1944, 1954 and 1959 for
female. Besides, the upper limit of life, that is Z , is estimated between 115 and 262, while those are estimated
between 117 and 165 for females with sufficiently small p-values.
Table 1. Estimates for male (left) and female (right)
Year
1948
1953
1958
1963
1968
1973
1978
1983
1988
1993
Gamma
-0.16
-0.13
-0.22
0.12
-0.08
-0.01
-0.07
-0.09
-0.13
-0.10
P-value
(one side)
0.08
0.13
0.01
0.88
0.09
0.43
0.04
0.00
0.00
0.00
P-value
(one side)
117.06
0.00
115.77
0.00
108.41
0.00
Omega
127.38
262.00
129.43
123.44
116.62
122.40
Year
1944
1949
1954
1959
1964
1969
1974
1979
1984
1989
0.00
0.00
0.00
0.00
0.00
0.00
Gamma
0.18
-0.07
0.07
0.03
-0.09
-0.12
-0.03
-0.05
-0.12
-0.12
P-value
(one side)
0.85
0.16
0.88
0.70
0.01
0.00
0.11
0.00
0.00
0.00
Omega
P-value
(one side)
128.71
0.00
121.78
117.83
165.00
143.40
120.00
120.83
0.00
0.00
0.00
0.00
0.00
0.00
References
Arssen, K., de Haan, L. (1994). On the maximal life span of humans. Mathematical Population Studies, Vol. 4,
259-281.
Kannisto, V. (1999). Trends in the mortality of the oldest-old, Statistics, Registries, and Science: Experiences from
Finland (ed. Alho, J.), Statistics Finland, 177-194.
Keiding, N. (1990). Statistical inference in the Lexis diagram. Philosophical Transactions of the Royal Society of
London Series A, Vol. 332, 487-509.
Wilmoth, J.R., Andreev, K., Jdanov, D. and Glei, D.A. (2007). Methods Protocol for the Human Mortality
Database. Vachonhttp://www.mortality.org/Public/Docs/MethodsProtocol.pdf.
Acknowledgment. This work was supported by Grant-in-Aid for Scientific Research 24500346.
大規模経済系データにおける様々な多重代入法アルゴリズムの検証
高橋
将宜
伊藤 孝之
はじめに1
データが欠測している場合、利用可能なデータサイズが縮小し、効率性が低下する。さら
に、観測値と欠測値との間に体系的な差異が存在する場合、統計分析の結果に偏りが発生す
るおそれがある。したがって、実際の統計分析においては、何らかの形で欠測値に対処する
ことがほとんど常に必須なことであり、欠測データの対処法として多重代入法(Multiple
Imputation)2が提唱されてきた(Rubin, 1987)。
多重代入法と一口に言っても、ソフトウェアに実装されているアルゴリズムには様々な方
法があり、マルコフ連鎖モンテカルロ法(MCMC: Markov chain Monte Carlo)、完全条件付指定
(FCS: Fully Conditional Specification)、EMB(Expectation-Maximization with Bootstrapping)の 3 つ
を有力なものとして挙げられる。現時点において、3 つのアルゴリズム間の優劣は判然とし
ていない。また、多重代入法の擬似データ数(M)をいくつに設定すればよいかについて、明確
な答えは見つかっていない。そこで、本稿では、大規模データセットとしての経済センサス
‐活動調査の速報データを用いて、3 つのアルゴリズムの比較検証を行った。さらに、経済
センサス‐活動調査の速報データに基づくシミュレーションデータを用い、多重代入法の擬
似データ数(M )をいくつに設定すればよいかについて検証を行った。
1.
多重代入法の理論
多重代入法では、観測データを条件として、欠測データの事後分布を構築し、この事後分
布からの無作為抽出を行うことで、補定にまつわる不確実性を反映させた M 個(M > 1)の補定
済データセットを生成することにより、欠測値を M 個のシミュレーション値に置き換える。
これら M 個の補定済データセットを別々に使用して統計分析を行い、しかるべき手法により
結果を統合し、点推定値を算出する。
2.
多重代入法の M 数
シミュレーションでは、一般的に、数百以上の副標本(M > 100)を生成する必要があり、コ
ンピュータの能力が許す限り多くの繰返しを行うべきだと考えられるが、元来、Rubin (1987)
によると、多重代入法の M は非常に小さい数字で十分だとされている。一方、近年では、M
数に関して Rubin (1987)への反論が展開されているものの、十分な結論を得るにいたっていな
い(Hershberger and Fisher, 2003; Carpenter and Kenward, 2007; Bodner, 2008) 。
独立行政法人統計センター統計情報・技術部統計技術研究課上級研究員
独立行政法人統計センター製表部管理企画課経済センサス業務推進室統計専門職
1
本研究の分析結果は、総務省・経済産業省『平成 24 年経済センサス‐活動調査』の速報結果の調査票情報を基
に著者が独自集計したものである。また、本稿の内容は、筆者の個人的見解を示すものであり、機関の見解を示
すものではない。
2 「多重代入法」とは、Multiple Imputation の訳である。総務省統計局及び統計センターでは、Imputation の訳語
として「補定」を用いているが、Multiple Imputation の訳語としては「多重代入法」の呼び名が一般的に流通し
ている(高橋, 伊藤, 2013, p.20)。よって、本稿においても、
「多重代入法」の用語を用いる。
1
3.
データセットと分析
3.1 経済センサス‐活動調査の速報データを用いた分析
2012 年 2 月に実施された経済センサス‐活動調査の速報データ（産業大分類 I の単独事業
所（個人経営以外）
）のデータ（観測数 277,263）を用い、R の多重代入法パッケージ Amelia
(EMB)、MICE (FCS)、Norm (MCMC)の比較検証を行った。欠測の発生メカニズムは、MAR
に基づき、売上高（自然対数）データの 20％（55,500 個）を人工的に欠測させた。また、資
本金（自然対数）データの 5％（13,600 個）を無作為に人工的に欠測させ（MCAR）、事業従
事者数（自然対数）には欠測を発生させていない（欠測率 0％）
。シードを 100 個生成し、3
つの多重代入法プログラムにより分析を行った。分析結果の詳細は、当日報告する。
3.2 シミュレーションデータを用いた多重代入データセット数 M の検証
自然対数変換した経済センサス‐活動調査の速報データの情報（平均値、分散・共分散な
ど）を基に、多変量正規分布によって観測数 1000、3 変量のシミュレーションデータセット
を 3 つ生成した。シミュレーションデータの基になったデータは、3.1 項と同じく、産業大分
類 I の単独事業所（個人経営以外）を用いた。Amelia、MICE、Norm において、シミュレーシ
ョンデータセット（欠測率 = 10%, 20%, 30%, 40%, 50%, 60%, 70%, 80%）に多重代入（M = 2, 5,
10, 15, 20, 30, 40, 50, 60, 70, 80, 90, 100, 200, 300, 400, 500）を施し、多重代入済データセットを
とその標準誤差の推定を行った。分析結
用いて、 ቀ売上高
ቁൌȽ
ෝ Ⱦ൫事業従事者数 ൯におけるȾ
୧
୧
果の詳細は、当日報告する。
4.
結語
本研究では、Amelia と MICE の補定値の精度には大きな差がないことが分かった。一方、
Amelia の計算処理速度は極めて速く、Norm は 27 万 3 変量のデータセットを分析すること
ができなかった。多重代入法の擬似データ数 M については、概ね 5∼10 では少なすぎ、20∼
50 程度が適切だと考えられる。欠測率に応じて、20%未満ならば M = 20、20%∼30%ならば
M = 30、30%∼40%ならば M = 40、40%∼50%ならば M = 50 といった具合に設定することが
適切だと思われる。また、欠測率に関わらず、M = 100 を超えて得られるものは非常に少ない。
欠測率が 50%を超え始めると、たとえ M 数を数百まで拡大したとしても、補定値の精度を保
証できなくなると考えられる。
参考文献
[1]
[2]
[3]
[4]
[5]
Bodner, Todd E. (2008). ³What Improves with Increased Missing Data Imputations?,´ Structural Equation Modeling
vol.15, pp.651-675.
Carpenter, James R. and Michael G. Kenward. (2007). Missing Data in Clinical Trials²A Practical Guide.
Birmingham: UK National Health Service, National Co-ordinating Centre for Research on Methodology.
Hershberger, Scott L. and Dennis G. Fisher. (2003). ³A Note on Determining the Number of Imputations for Missing
Data,´ Structural Equation Modeling vol.10, no.4, pp.648-650.
Rubin, Donald B. (1987). Multiple Imputation for Nonresponse in Surveys. New York: John Wiley & Sons.
高橋将宜, 伊藤孝之. (2013). 「経済調査における売上高の欠測値補定方法について∼多重代入法による精度
の評価∼」, 『統計研究彙報』第 70 号 no.2, 総務省統計研修所, pp.19-86.
2
国勢調査ミクロデータを用いた匿名化技法の有効性の検証*
明海大学経済学部 伊藤 伸介**
(独)統計センター 星野 なおみ
1．はじめに
わが国においては、平成 25 年末に平成 12 年国勢調査の匿名データの作成・提供が開始され
た。その一方で、将来的には、小地域分析用の匿名データ等、別のタイプの国勢調査の匿名デー
タの要望が出てくる可能性があり、その予備的な研究としてミクロデータに対する匿名化技法の
適用可能性を検証することは有用であると考えられる。そこで、本稿では、各種匿名化技法を用
いて作成された国勢調査の匿名化ミクロデータを対象に、秘匿性と有用性に関する実証研究を行
うことによって、匿名化技法の有効性の検証を行った。
2. 本研究における匿名化ミクロデータの作成方法
本研究では、国勢調査を例に匿名化ミクロデータを試行的に作成した。本研究で使用するデー
タは、平成 17 年国勢調査(以下「国調」と略称)の個票データをもとに特定の地域(以下「地域 A」
と呼称)のレコードから作成したテストデータ(約 100,000 レコード)である。このテストデータ
には、個人単位で抽出した一般世帯の世帯主のレコードのみが含まれている。
匿名化ミクロデータの作成手順は以下のとおりである。(1)国調のテストデータに対して、労
働力状態、従業上の地位と年齢についてはリコーディングを、世帯人員と年齢に関してはトップ
コーディングをそれぞれ適用した。(2)リコーディングとトップコーディングを施したデータに
対して様々な標本抽出率によるサンプリングを行った(本稿では、サンプリング率が 10%の場合
の結果に基づいて議論する)。(3)(1)と(2)の匿名化技法が適用されたデータに対して、スワッピン
グを適用した。具体的なスワッピングの手順としては、(1)世帯主との続き柄等の 11 のキー
変数を用いて標本一意を計測した上で、スワッピングの対象となるレコードを選出し、(2)
スワッピングの対象レコードの中で優先度の高いレコードをスコアに基づいて探索し、(3)
対象レコードに対してターゲット・スワッピングとランダム・スワッピングを適用している。
なお、スワッピング率については、1%、2%、3%、5%、10%、20%、30％が設定されている。
3. 国勢調査ミクロデータを用いた秘匿性の評価―マッチングの試み―
本研究では、秘匿性の評価方法の１つとして、外部情報とミクロデータのマッチングに焦点を
当てる。具体的には、国調以外の政府統計のミクロデータを外部情報とみなした上で、匿名化ミ
クロデータと外部情報とのマッチングを試みた。本研究においては、スワッピングが施された
10％抽出の国調の匿名化ミクロデータに対して、平成 20 年住宅・土地統計調査 (以下「住調」
と略称)の地域 A に該当するレコードを含む個票データ(約 10,000 レコード)とのマッチングの実
験を行った。住調については、国調のテストデータと調査区が重複するように対象レコードが選
定される。
本実験では、国調と住調の両方の共通の調査事項から選ばれた県市区町村番号、世帯人員、性
別といった変数を含む様々なキー変数の組み合わせに基づいて、国調と住調のマッチングを行っ
た。本実験においては、国調と住調においてキー変数の区分を合わせた上で、マッチングを行う。
また、本実験では、国調の匿名化ミクロデータにおいて上記のキー変数で一意になったレコード
を対象に、住調とのマッチングを行っている。さらに、マッチングにおける国調と住調の調査年
次の違いについては、国調のレコードにおいて年齢を加算することによって、年次の調整がなさ
れている。
国調の匿名化ミクロデータと住調の個票データのマッチングに関する実験結果から、主として、
以下の点が明らかになった。第 1 に、キー変数の数が多くなるにつれて、国調の匿名化ミクロ
データにおける標本一意の数が増大していることがわかるが、この傾向は、ターゲット・スワッ
ピングだけでなく、ランダム・スワッピングの場合にも当てはまることが確認された。また、タ
ーゲット・スワッピングとランダム・スワッピングのいずれにおいても、スワッピング率が高く
なるにつれて、国調の匿名化ミクロデータにおける標本一意の数が小さくなっている。その理由
*
本稿の内容は個人的な見解を示すものであり、統計センターの見解を表すものではないことに留意され
たい。
** (独)統計センター非常勤研究員
1
としては、スワッピングによって、標本一意に該当するレコードが度数 2 以上のセルに該当す
るグループに移動したことが考えられる。
第 2 に、スワッピング率が上昇するにつれて、マッチングされたレコードの中でスワッピン
グされたレコードの比率が高くなっているものの、スワッピング率が 30％の場合でもその比率
は最大で約 60％となっている。このことは、本実験においてスワッピング率を上げても、国調
の匿名化ミクロデータのレコードの中で住調の個票データとマッチングされたレコードの一部
に対してのみ、スワッピングが施されていることを意味する。
第 3 に、本実験では国調の標本一意に対する「真のマッチング」1の比率が算出されているが、
標本一意に占める真のマッチングの比率を確認すると、ターゲット・スワッピングおよびランダ
ム・スワッピングのいずれについてもその比率は約 1％∼2％あって、非常に低くなっている。
このことは、国調と住調においてマッチングされたレコードの組については、同じ調査区番号を
有していてもそれらが同一の世帯のレコードに該当する可能性が低いことを示している。
4. スワッピングの有効性の評価に関する研究
本節では、主として「特殊な一意(special uniques)」のレコードを対象に適用されたスワッピ
ングの有効性を検証するために、スワッピングが適用された匿名化ミクロデータ(以下「スワッ
ピング済データ」)を対象に、その有用性と秘匿性の評価を行った。本研究では、Shlomo et
al.(2010)に基づいて，有用性と秘匿性の評価指標を作成した。有用性の評価指標(DU)は、絶対
距離の平均値(average absolute distance)で与えられる。具体的には、原データとスワッピング
済データの両方で集計値を作成した上でセルごとの度数の差の絶対値に関する平均値が算出さ
れる。一方、秘匿性の評価指標(DR)については、原データにおけるクロス表の中で度数 1
であるセル数に対して、スワッピング済データにおけるクロス表の中で度数 1 でありかつ
スワッピングされていないセル数の比率が算出されている。
これらの有用性と秘匿性の評価指標に基づいて、本研究は、R-U マップ(R-U Confidentiality
Map)を作成し、ターゲット・スワッピングとランダム・スワッピングについて、有用性と秘匿
性の相対比較を試みた。なお、R-U マップで使用する有用性と秘匿性の評価指標に関しては、
キー変数の中のあらゆる 2 変数の組み合わせについて計算された評価指標の平均値がそれぞれ
用いられている。
R-U マップの結果から主として以下の点が明らかになった。第 1 に、スワッピング率が上昇
するほど、情報量損失値としての DU の数値が大きくなり、DR の数値が小さくなる傾向にある
ことが確認できる。第 2 に、ランダム・スワッピングと比較して、ターゲット・スワッピング
における DU の数値が大きく、DR の数値が小さくなることが確認できる。このことは、ターゲ
ット・スワッピングについては、ランダム・スワッピングと比較した場合、有用性は低いものの
秘匿性は相対的に高いことを示唆している。
5. むすびにかえて
本稿では、国勢調査のミクロデータを用いて、スワッピング等の各種匿名化技法の有効性の検
証を行った。外部情報とのマッチングの試みとして、国調の匿名化ミクロデータと住調の個票デ
ータとのマッチングを行ったが、マッチングの精度は高くないことから、外部情報とのマッチン
グの可能性という観点から見た場合、露見のリスクは低いとみなすことができる。また、本研究
では、R-U マップをもとに、スワッピング済データにおける有用性と秘匿性の検証も行った。
有用性あるいは秘匿性に関する閾値を設定することができれば、適切なスワッピング率およびス
ワッピングの方法を選択することが可能になることがわかった。
参考文献
Shlomo, N., Tudor, C., Groom, P. (2010)
Data Swapping for Protecting Census Tables ,
Domingo-Ferrer, J. and Magkos, E.(eds) Privacy in Statistical Databases UNESCO Chair in Data
Privacy International Conference, PSD 2010 Corfu, Greece, September, 2010 Proceedings, Springer,
pp.41-51.
1本研究における真のマッチングとは、国調の匿名化ミクロデータと住調の個票データの間で
1 対 1 にマッ
チングされたレコードにおいて、調査区番号も同一であることが確認されることである。マッチングされ
たレコードが真のマッチングに該当するレコードかどうかを確認するために、本実験では、統計センター
の内部資料に基づいて、国調と住調でマッチングされたレコードの組を対象に、それらのレコードの調査
区番号が同一か否かについて検証作業を行った。
2
ピットマン確率分割に付随する極限分布への一考察
鹿児島大学 大和 元
1 始めに
Pitman’s sampling formula に従う確率分割について、分割の要素の個数 Kn の漸近分布に関連
した事項を考えるが、Kn の分布は次で与えられる。
P (Kn = k) =
θ[k:α]
c(n, k, α)α−k
θ[n]
ここで、c(n, k, α) = (−1)n−k C(n, k, α) であり、 C(n, k, α) は C-number と言われ次式の右辺の
係数である：(st)(n) =
!n
k=1
C(n, k, s)t(k) 、ただし、t(k) = t(t − 1) · · · (t − k + 1)。
補題 1.1 (Pitman (1996), Yamato and Sibuya (2000))
µ"r = lim E
n→∞
"# K $r %
n
nα
#
# θ $[r] Γ(θ)
θ $[r] Γ(θ + 1)
= 1+
=
.
α
Γ(θ + 1 + rα)
α
Γ(θ + rα)
この µ"r は、密度関数が
Γ(θ + 1) θ
x α gα (x) (0 < α < 1, θ > −α)
Γ( αθ + 1)
(#)
で与えられる分布の r 次の積率である。ここで、gα はパラメータ α の Mittag-Leffler distribution
(M L(α)) の密度関数である。 密度関数 gα (0 < α < 1) は次式を満たす一意な関数である。
& ∞
Γ(p + 1)
xp gα (x)dx =
(∀ p > −1)
Γ(pα + 1)
0
(Pitman (1996) p.5 and Pitman (2002), p.20-21).
密度関数 (#) を持つ分布を Pitman-Mittag-Leffler distribution と言う事にし、P M L(θ, α) で表
す。 P M L(0, α) = M L(α) である。法則収束を → で表すと、補題 1.1 から次が得られる。
d
補題 1.2 (Pitman (1996), Pitman (2002))
Kn d
→ P M L(θ, α)
nα
as
n→∞
分布 P M L(θ, α) を、従って本質的には M L(α) を、より分かり易い形で表すことが本研究の目
的である。以下に示す様に、パラメータを特別な場合に限れば表す事が出来た。
1
2 P M L の分布
命題 2.1 |X| (X ∼ N (0, 2)) は Mittag-Leffler 分布 M L(1/2) に従い、
E(|X|p ) =
Γ(p + 1)
Γ( 12 p + 1)
(p > −1).
命題 2.2 r = 1, 2, · · · について、X1 ,· · · , Xr は独立で、N (0, 2) に従うとする。
1
1
1
Yr = |X1 | · |X2 | 2 · |X3 | 22 · · · |Xr | 2r−1
は Mittag-Leffler 分布 M L(1/2r ) に従う。 或は、 r = 1, 2, · · · について、Z1 ,· · · , Zr は独立で、
N (0, 1) に従うとする。
1
r
1
1
Yr = 21−1/2 |Z1 | · |Z2 | 2 · |Z3 | 22 · · · |Zr | 2r−1
は Mittag-Leffler 分布 M L(1/2r ) に従う。
Y2 = |X1 | · |X2 | 2 の密度関数は
1
2
f|X1 |·|X2 |1/2 (x) =
π
&
∞
0
exp
'
−
$(
1 # x2
4
+
t
dt (x ≥ 0)
4 t2
自由度 ν のカイ２乗分布の密度関数と分布関数を hν (x)、Hν (x) とおく。
命題 2.3 パラメータ α = 1/2、 θ > −1/2 について、Pitman-Mittag-Leffler 分布 P M L(θ, α) の
分布関数と密度関数は、
Fθ,1/2 (x) = H2θ+1
)
*
x2
,
2
fθ,1/2 (x) = xh2θ+1
)
*
x2
,
2
x > 0.
θ > −α, α = 1/4 に対して、P M L の密度関数は、
&
' 1 # x2
$(
Γ(θ + 1) θ 2 ∞
4
α
fθ,α (x) =
x
exp
−
+
t
dt
π 0
4 t2
Γ( αθ + 1)
α '= 1/2r (r = 1, 2, · · · ) の場合の P M L の分布を調べることが、今後の問題である。
参考文献
Pitman, J. (1996), Notes on the two parameter generalization of Ewens’ random partition structure (unpublished manuscript).
Pitman, J. (2002), Poisson-Kingman Partitions, Technical Report 625, Dept. Statistics,
U.C. Berkeley.
Yamato H. and Sibuya, M. (2000), Moments of some statistics of Pitman sampling
formula, Bull. Inform. Cybern., 32, 1–10
2
:
2013-11-27/29
S ∈ Pn,k |(S ∈ Pν,κ , n < ν)
(ν, κ)
.
, S ∈ Pn,k
. n :=
Kn :=
!k
ci , 1 ≤ i ≤ k,
.
{c1 , . . . , ck }
Ci
i=1 ci ,
k,
,
EPSF (n; θ, α)
(ci = j)
Ci
!n
j=1 Sj , (Sj
,
.
)
p.m.f., P{Kn = k} =: fn (k)
n − kα
θ + (k − 1)α
fn (k) +
fn (k − 1), 1 ≤ k ≤ n,
θ+n
θ+n
1
fn (k) =
Sn,k (θ| − 1)n (θ| − α)k , 1 ≤ k ≤ n,
(θ| − 1)n
fn+1 (k) =
Sn,k := Sn,k (−1, −α, 0)
Sn,k (a, b, c)
(n; θ, α)
.
EPSF (θ, α)
(1)
(2)
EPSF-K
.
#
$
n
"
1
1 (1 − α| − 1)j−1 sj
w(s; n, k) := P{S = s|(S ∈ Pn,k )} =
.
Sn,k (−1, −α, 0)
sj !
j!
(3)
j=1
.
θ
,
.
.
.
Kn
fn (k)
.
,
{(n, k); 1 ≤ k ≤ n < ∞}
. Kν
,
(3)
Kν−1
P{Kν−1 = κ − 1|Kν = κ} = Sν−1,κ−1 /Sν,κ
P{Kν−1 = κ|Kν = κ} = (ν − 1 − κα)Sν−1,κ /Sν,κ = 1 − Sν−1,κ−1 /Sν,κ .
n
P{Kn = k}
.
gn (k) = (n − kα)
Kν
gn (k)
(
)
.
Sn,k
Sn,k
gn+1 (k) +
gn+1 (k + 1),
Sn+1,k
Sn+1,k+1
1 ≤ k ≤ n.
S|(S ∈ Pν,κ )
p.m.f., gn (k), 1 ≤ k ≤ n
gk (n) = gk (n; ν, κ, α)
gκ (ν; ν, κ, α) = 1
g1 (1; ν, κ, α) = 1
%[ν, κ] := {(n, k); max(1, n − ν + κ) ≤ k ≤ min(κ, n)}
(4)
,
gk (n; ν, κ, α) = 0, k > min(κ, n),
gn (k) :=
&
Sn,1 /Sn+1,1 = 1/(n − α).
(4)
.
,
,
n
1
E(Kn |S ∈ Pν,κ ) =
Sν,κ
(4)
:
min(κ,ν−n)
%
k=1
## $ #
$$
ν
ν−n
Sν−k,κ−1 (1−α|−1)k−1
−
,
k
k
(µ, λ) ∈ ![ν, κ]
(ν, κ)
.
,
, t
,
.
κ
. α
,
ν = n0 + 1, n0 + 2, . . .
k0 ≤ κ ≤ k0 + ν − n0
1 < n < ν. (5)
t
.
t ∈ Pn0 ,k0
S|(S ∈ Pν,κ )
, ν
. (n0 , k0 ) ∈ %[ν, κ]
κ̂ = argmaxλ gn0 (k0 ; ν, λ)
κ̂
k0
,α
(b)
.
MLE
.
Osa Peninsula, Costa Rica,
2
. (Colwell, et al. 2012). (a) (n, k) = (976, 140) , (b) (n, k) = (237, 112) ,
(a)
237
69
.
(b)
n = 976
gn (k)
234
.
.
[1] Colwell, R.K. et al.(2012) Models and estimators linking individual-based and sample-based rarefaction,
extrapolation and comparison of assemblages, Journal of Plant Ecology The Annals of Applied Probability ,
15-1, 3–21.
[2] Gnedin, A. and Pitman, J. (2006) Exchangeable Gibbs partitions and Stirling triangles, Mathematical Sciences, 138, 5674–5685. (original —Russian version: Zapiski Nauchnnykh Seminarov ROMI, 325, 2005,
83–102.)
[3] Lijoi, A., Prünster, I. and Walker, S.G. (2008) Baysian nonparametric estimators derived from conditional
Gibbs structures, The Annals of Applied Probability, 18-4, 1519–1547.
[4] Pitman, J. (2006) Combinatorial Stochastic Processes, Lecture Notes in Mathematics, 1875, Springer, New
York, NY.
[5] Sibuya, M. (2013) Prediction in Ewens-Pitman Sampling Formula and random samples from numberpartitions, Annals of the Institute of Statistical Mathematics, (in print).
Semiparametric Efficiency for the
Quantile-Regression-Based L-Estimation
Takayuki Shiohama (Tokyo Univerisity of Science)
Hiroyuki Taniai (Waseda Univerisity)
Abstract
We study a semiparametrically efficient estimation for a quantile-regressionbased L-estiamtors. The L-estimators are based on linear combinations of the linear
quantile regression estimator, introduced by Koenker and Bassett (1978), and our
proposed L-estimators are based on the one-step estimator of the type discussed by
Hallin and Werker (2003) and Hallin et al. (2008). In constructing the estimator,
we use the residual signs and ranks based optimal score function with a general
reference density. The asymptotic properties of the proposed estimators are discussed. Their finite sample properties are illustrated through a set of Monte Carlo
simulations and an empirical applications.
1.
Introduction
Quantile regression introduced by Koenker and Bassett (1978) is similar to the standard
linear regression, expect that we replace conditional expectations E(Y |X) by a conditional quantile FY−1|X (τ ). Quantile regression provides us a more complete summary of the
conditional distribution. A throughout review of QR can be found in Koenker (2005).
Linear combinations of order statistics, called L-estimators, have long been of interest as estimates for the classical location model for robust alternatives. The problem of
constructing asymptotically efficient estimators of a finite dimensional Euclidian parameter θ in a semi-parametric model with an infinite dimensional nuisance parameter F has
attracted considerable attention throughout last two decades.
2.
Quantile-Regression-based L-Estimation
Let FY |X (y) and FY−1|X (τ ) = inf{y : FY |X (y) ≥ τ } denote the conditional distribution
function and the τ -quantile of Y given X, respectively. The model of our concern is
represented as observations {yt }nt=1 from statistical experiments (Y n , B(Y n ), PYn ) satisfying
Yt = Xt" β(τ ) + ε(τ )t , Fε−1 (τ ) = 0,
with some covariate vector Xt = [1, X2,t , . . . , Xd,t ]" , for t = 1, . . . , n, and an unknown
regression parameter to be estimated β ∈ Rd . Equivalently, we have by taking the τ quantile on the both hand side of (2.1), conditionally onXt , we obtain the value of β in
the relation FY−1|X (τ ) = Xt" β(τ ). This is known as basic conditional quantile regression
model.
In this paper, we consider the asymptotic properties of L-estimation of the linear
models via quantile regression and its one-step estomator. Our discussion focuses on the
semiparametrically efficient estimation which is based on the quantile regression process
β!nQR (·) dedined, following Koenker and Bassett (1978), as
β!n(QR) (τ ) = argmin
b∈Rd
n
"
t=1
ρτ (Yi − Xi" β),
1
τ ∈T,
where ρτ (u) = u(τ − I(u < 0)) is called check function.
Let J : (0, 1) → R be the weight function. For the parameter of interest is
# 1
$# 1
%
"
"
φJ ◦ ψτ (PY ) :=
J(τ )X β(τ )dτ = X
J(τ )β(τ )dτ ,
0
and its QR based estimator is, for τi =
0
i
,i
n+1
= 1, . . . , n,
n
1"
φJ ◦ ψτ (Pn,Y ) :=
J (τi ) X " β!n(QR) (τi )
n i=1
The procedure that we will apply here to achieve semiparametric efficiency is based on
the invariance principle, as introduced by Hallin and Werker (2003). By following Hallin
and Werker (2003), a semiparametrically efficient procedure can be obtained by projecting
(n)
(n)
∆bn ;f on some σ-field to which the generating group for {Pbn ,f |f ∈ F α } becomes maximal
invariant. For the quantile-restricted regression model, such a σ-field is studied by Hallin
et al. (2008) and found to be generated by signs and ranks of the residuals.
(n)
Here, let us denote the sign of a residual as Sbn ,i , the rank of a residual as Rbn ,i , and
the σ-field generated by them as
(n)
(n)
(n)
Bbn := σ(Sbn ,1 , . . . , Sbn ,n ; Rbn ,1 , . . . , Rbn ,n ).
A “good” inference should be based on
n
1 " (n)
(n)
(n)
(n)
(n)
(n)
∆β,f := Eβ,f [∆β,f |Bβ ] = 1/2
Eβ,f [ϕf (F −1 (Uβ,i ))|Bβ ]Xi
n
&
i=1
where Uβ,i := F (ξβ,i ) = F (Yi − β T Xi ) is the i.i.d. uniform on [0, 1] under Pβ,f .
Then we propose a one-step estimator for the quantile regression coefficient processes
and discuss its semiparametric inference. From Theorem18.6 of Kosorok (2008), we need
to show that the composite function φ is Hadamard differentiable, which yields our proposed one-step L-estimators are semiparametrically efficient.
References
Hallin, M., Vermandele, C., and Werker, B. J. M. (2008). Semiparametrically efficient
inference based on signs and ranks for median-restricted models. J. R. Stat. Soc. Ser.
B Stat. Methodol., 70(2):389–412.
Hallin, M. and Werker, B. J. M. (2003). Semi-parametric efficiency, distribution-freeness
and invariance. Bernoulli, 9(1):137–165.
Koenker, R. (2005). Quantile regression, volume 38 of Econometric Society Monographs.
Cambridge University Press, Cambridge.
Koenker, R. and Bassett, Jr., G. (1978). Regression quantiles. Econometrica, 46(1):33–50.
Kosorok, M. R. (2008). Introduction to empirical processes and semiparametric inference,
Springer Series in Statistics. Springer, New York.
2
On a tangent space for the coefficient functions
of Quantile Regression.
谷合 弘行 (Hiroyuki TANIAI, 早稲田大学)
シンポジウム「統計科学の新展開」 2013.11/27-29, 金沢大学サテライト・プラザ
本発表の概要: 一般の分位数回帰モデルに従う確率変数は、一様確率変数を引数とするランダム係
数（random coefficient）表現を持つ。本発表では、これを一種の混合分布モデル（mixture model）
とみなすことによって無限次元パラメータである係数関数の接空間を導くことを試みた。また、さ
らにこの結果を用いた漸近有効な推定量の構成についても議論した。以下がその概要である：
We observe 𝑿𝑡 = (𝑌𝑡 , 𝒁𝑡 ), 1 ≤ 𝑡 ≤ 𝑛 where 𝑌𝑡 is a random variable and 𝒁𝑡 is a 𝑑-dimensional
random vector, satisfying 𝐹𝑌−1 (𝜏 ∣𝒁𝑡 = 𝒛𝑡 ) = 𝒛𝑡⊤ 𝜷(𝜏 ), for all 𝜏 ∈ (0, 1). Here, we denote the
conditional quantile by 𝐹𝑌−1 (𝜏 ∣𝑆) := inf{𝑦 : 𝑃 𝑌 {𝑌 ≤ 𝑦∣𝑆} ≥ 𝜏 } with some set 𝑆 given. This is
the QR model (cf. Koenker (2005) and/or Koenker and Xiao (2006)) so has an interpretation
as the following Random Coefficient Regression model:
d
𝑌𝑡 = 𝒁𝑡⊤ 𝜷(𝑈 ),
𝑈 ∼ Uniform[0,1],
(1)
where 𝜷(⋅) : (0, 1) → ℝ𝑑 is the vector of non-decreasing functions. In particular, if we assume
the 1st entry of 𝒁𝑡 is 1 and set 𝜷(𝑢) = (𝜃1 + 𝑏(𝑢), 𝜃2 , . . . , 𝜃𝑑 )⊤ , then 𝜷(𝑢) is studied by the
location model with restriction to the innovation’s 𝑢-quantile to be 0. But the same 𝜷(𝑢) can be
studied in a scale model with 𝑢-quantile being 1, by setting 𝜷(𝑢) = (𝜃1 𝑐(𝑢), . . . , 𝜃𝑑 𝑐(𝑢))⊤ . So,
in order to have the efficiency gain, different perspective yields different regularity conditions.
In this talk, we focus on the modelling (1), where 𝑌 is generated by deterministic function 𝜷,
uniform r.v. 𝑈 , and “random nuisance parameter” 𝒁 (cf. Pfanzagl (1982)). The distribution
of 𝑿 is now be seen from
∫
𝑃𝜷,Γ (𝐴 × 𝐵) =
𝑃 𝑌 (𝐴; 𝜷∣𝒛)Γ(𝑑𝒛), 𝐴 ⊂ ℝ, 𝐵 ⊂ ℝ𝑑
𝐵
𝒁
where Γ := 𝑃 plays a role of mixing distribution. The crucial difference is that the parameter
of our interest, 𝜷(⋅), is infinite-dimensional.
To determine the tangent spaces for 𝜷 and for Γ, we define the subfamilies
𝔓Γ := {𝑃𝜷,Γ : 𝜷 ∈ ℳ𝑑 },
𝔔𝜷 := {𝑃𝜷,Γ : Γ ∈ G },
with ℳ being the set of non-decreasing functions on (0, 1) and G being the set of distributions
which 𝒁 might have. A tangent cone for 𝜷(⋅) (at 𝑃 ) is denoted by 𝑇𝑃 (𝔓Γ ), which is the set
of (score) functions 𝑔 satisfying 𝑃 [𝑔 2 ] < ∞ and
𝑑𝑃𝛿
= 1 + 𝛿(𝑔 + 𝑟𝛿 ),
𝑑𝑃
1
𝑃 [𝑟𝛿2 ] = 𝑜(𝛿 0 )
for the path 𝑃𝛿 → 𝑃 (in 𝔓Γ ) as 𝛿 → 0. Now, denote the true value of parameters (𝜷0 (⋅), Γ0 )
and think of a contamination 𝜷(⋅, 𝛿) such that 𝜷𝛿 (𝑢) := 𝜷(𝑢, 𝛿) → 𝜷0 (𝑢) as 𝛿 → 0. For given
𝑢 ∈ (0, 1), this 𝜷(𝑢, 𝛿) may differ at 𝛿 ∕= 0 from one parametrization to another. Also, as
∂
for the partial derivative ∂𝛿 𝜷(𝑢′ , 𝛿 ′ ) := ∂𝛿
𝜷(𝑢, 𝛿)∣(𝑢,𝛿)=(𝑢′ ,𝛿′ ) , we may think of the function
−1
0
∂𝛿 𝜷(⋅, 0) as 𝛿 (𝜷𝛿 − 𝜷0 )(⋅) + 𝑜(𝛿 ) so that each ∂𝛿 𝛽𝑖 (⋅, 0) ∈ 𝐿2 (𝑃𝜷0 ,Γ0 ) if 𝜷𝛿 is such that
√
𝑛-consistent. In fact, each ∂𝛿 𝛽𝑖 (⋅, 0) can be assumed square-integrable under 𝑃 𝑌 (⋅; 𝜷0 ∣𝒛) as
well.
Proposition 1 The tangent space for 𝜷(⋅) is found to be
𝑇𝑃𝜷0 ,Γ (𝔓Γ ) = {𝐴Γ 𝒉 : ℎ𝑖 ∈ 𝐿2 (𝑃 𝑌 (⋅; 𝜷0 ∣𝒛)), 𝑖 = 1, . . . , 𝑑}, where
𝑑 ( 𝒛 ⊤ 𝒉(𝑢) )$$
𝑑
(𝐴Γ 𝒉)(𝑦) := −
. (𝜷0′ (𝑢) =
𝜷0 (𝑢))
$
𝑑𝑢 𝒛 ⊤ 𝜷0′ (𝑢) 𝑢=(𝒛⊤ 𝜷0 )−1 (𝑦)
𝑑𝑢
On the other hand, those for Γ will be 𝑇𝑃𝜷,Γ0 (𝔔𝜷 ) = {(𝑦, 𝒛) 1→ 𝑘(𝒛) : ∃ 𝑘 ∈ 𝑇Γ0 (G )}.
This 𝐴Γ : 𝐿2 (𝑃 𝑌 (⋅; 𝜷0 ∣𝒛))𝑑 → 𝐿2 (𝑃𝜷0 ,Γ ) is the score operator for 𝜷(⋅). So the score operator
for the full model seems to be 𝐴(𝒂, 𝒃) ≡ 𝐴Γ 𝒂 + 𝒃⊤ ∇𝑘(𝒛), for (𝒂, 𝒃) ∈ 𝐻 := 𝐿2 (𝑃 𝑌 (⋅; 𝜷0 ∣𝒛))𝑑 ×
ℝ𝑑 . Then its adjoint 𝐴∗ : 𝐿2 (𝑃𝜷0 ,Γ0 ) → 𝐻 can be formally given by
(
[
)
)]/ 𝑌
𝒛
∂ ( ⊤ ′
𝐴∗ 𝑓 = 𝑃𝜷0 ,Γ0 ⊤ ′
𝒛 𝜷0 (𝑢0 ) ⋅ 𝑓 (𝒙)
𝑑𝑃 (𝑦; 𝜷0 ∣𝒛), 𝑃𝜷0 ,Γ0 [∇𝑘(𝒛)𝑓 (𝒙)]
𝒛 𝜷0 (𝑢0 ) ∂𝑦
with 𝑢0 = (𝒛 ⊤ 𝜷0 )−1 (𝑦). Using these operators, the efficient influence function for each 𝛽𝑖 (𝑢),
˜ 𝑖⋅ (𝑦) as in eq.(25.30) of van der Vaart (1998).
𝑖 = 1, . . . , 𝑑, can be obtained as 𝐴(𝐴∗ 𝐴)−1 𝝌
˜ 𝑖⋅ (𝑦) = (Γ[∂𝛿 𝜷(𝑢, 0)∂𝛿 𝜷(𝑢, 0)⊤ /𝒛 ⊤ 𝜷0′ (𝑢)]−1 Γ[∂𝛿 𝛽𝑖′ (𝑢, 0)∂𝛿 𝜷(𝑢, 0)], 0𝑑 ). More
In our case, 𝝌
specifically, we have the following result:
Proposition 2 The efficient function for each 𝜷𝑖 (𝑢) is found to be 𝜓˜𝑖 = 𝐴Γ 𝒌𝑖⋅ where 𝒌𝑖⋅ (𝑢) =
(𝑘𝑖1 (𝑢), . . . , 𝑘𝑖𝑑 (𝑢))⊤ is the solution of the following ODE
{ (𝒛 ⊤ 𝜷 ′′ (𝑢))2 − 𝒛 ⊤ 𝜷 ′′′ (𝑢)
𝒛 ⊤ 𝜷 ′′ (𝑢) ′
0
0
′′
𝑘𝑖𝑗
− ⊤ 0′
𝑘𝑖𝑗 +
}𝑘𝑖𝑗
𝒛 𝜷0 (𝑢)
𝒛 ⊤ 𝜷0′ (𝑢)
[ ∂ 𝜷(𝑢, 0)∂ 𝜷(𝑢, 0)⊤ ]−1
𝛿
𝛿
= 𝜷0′⊤ (𝑢)𝜷0′ (𝑢)∂𝑢𝛿 𝛽𝑗 (𝑢, 0)Γ
∂𝛿 𝛽𝑗 (𝑢, 0)
𝒛 ⊤ 𝜷0′ (𝑢)
with conditions 𝑘𝑖𝑗 (0) = 𝑘𝑖𝑗 (1) = 0.
As for the construction of the estimator which realizes this efficient influence function, we
made some conjectures and left it as a future work.
参考文献
Koenker, R. (2005). Quantile regression, volume 38 of Econometric Society Monographs.
Cambridge University Press, Cambridge.
Koenker, R. and Xiao, Z. (2006). Quantile autoregression. J. Amer. Statist. Assoc.,
101(475):980–990.
Pfanzagl, J. (1982). Contributions to a general asymptotic statistical theory, volume 13 of Lecture Notes in Statistics. Springer-Verlag, New York. With the assistance of W. Wefelmeyer.
van der Vaart, A. W. (1998). Asymptotic Statistics. Cambridge Series in Statistical and
Probabilistic Mathematics. Cambridge University Press, Cambridge.
2
2
t
1
1.
R×R
pij
, (i, j)
(i = 1, . . . , R; j = 1, . . . , R).
(Bowker, 1948)
:
pij = pji
(i < j).
.
(Agresti, 1983)
:
pij = θj−i pji
.
(i < j).
j −i
θ=1
.
2
U
V
,
E(U ) = µ1 , E(V ) = µ2 ,
2
Var(U ) = Var(V ) = σ , Corr(U, V ) = ρ.
:
f (u, v)
,
f (u, v)
= δ v−u ,
f (v, u)
,
!
"
µ2 − µ1
δ = exp 2
.
σ (1 − ρ)
,
2
,
,
.
Yamamoto et al. (2007)
.
,
2
1
t
.
278-8510
2641
e-mail: [email protected]
2. 2
t
m
2
t
U
V
,
E(U ) = µ1 ,
2
E(V ) = µ2 , Var(U ) = Var(V ) = mσ /(m − 2) (m > 2), Corr(U, V ) = ρ.
h(u, v)
,
$
%− m+2
2
1
Q(u, v)
#
h(u, v) =
1+
,
m
2πσ 2 1 − ρ2
,
&
'
1
2
2
Q(u, v) = 2
(u − µ1 ) − 2ρ(u − µ1 )(v − µ2 ) + (v − µ2 ) .
σ (1 − ρ2 )
h(u, v)
(
,
h(u, v)
2
)− m+2
:
2
(
)− m+2
− h(v, u)
= αm (v − u) (u < v),
2
) m+2
#
2(µ1 − µ2 ) (
2
2
αm =
2πσ 1 − ρ
.
mσ 2 (1 − ρ)
R×R
,
−
2
.
m (m > 2)
2
−
pij m+2 − pji m+2 = βm (j − i) (i < j),
βm
.
m
t
,
t
TS(m)
βm = 0
pji
.
.
.
−2/(m + 2)
.
pij
j−i
−2/(m + 2)
Agresti, A. (1983). A simple diagonals-parameter symmetry and quasi-symmetry model.
Statistics and Probability Letters, 1, 313-316.
Bowker, A. H. (1948). A test for symmetry in contingency tables. Journal of the American Statistical Association, 43, 572-574.
Yamamoto, H., Iwashita, T. and Tomizawa, S. (2007). Decomposition of symmetry into
ordinal quasi-symmetry and marginal equimoment for multi-way tables. Austrian
Journal of Statistics 36, 291-306.
順序カテゴリ正方分割表における対称性のモデルと分解
および併合した表に基づく対称性に関する尺度
東京理科大学大学院・理工学研究科 島田 文香
東京理科大学大学院・理工学研究科 田中 弥生
大阪大学・医学部
山本 紘司
東京理科大学・理工学部
富澤 貞男
第 1 部 Sum-symmetry モデルの提案とその分解
第一部では, モデルを 2 つ提案した。まず, 新たな対称構造をみるモデルとして
次のモデルを提案した: !
pij
(i,j)∈R(x)
=
!
(i,j)∈R(x)
pji (3 ≤ x ≤ 2r − 1),
ここに,
R(x) = {(i, j) | i + j = x, i < j},
である. ここではこのモデルを Sum-symmetry(SS) モデルと呼ぶ. 次に下記のモデ
ルを導入した:
!
pij
= ∆
(i,j)∈R(x)
!
(i,j)∈R(x)
pji (3 ≤ x ≤ 2r − 1).
このモデルを条件付き Sum-symmetry(CSS) モデルと呼ぶ. また, ∆ = 1 とおいた
CSS モデルは SS モデルである. このとき次の分解定理を得る:
定理 1: SS モデルが成り立つための必要十分条件は, CSS モデルと GS モデルの
両方が成り立つことである.
r × r 正方分割表において, (i, j) セル観測度数を nij とし, mij を対応するセル期待
度数とする. {nij } は多項分布に従うと仮定する. このとき, モデル M の適合度を
検定するための尤度比カイ二乗統計量は
2
G (M ) = 2
r !
r
!
nij log
i=1 j=1
"
nij
m̂ij
#
,
で与えられる. m̂ij はモデル M の下での mij の最尤推定値である. このとき, 次の
定理を得る:
定理 2. 検定統計量 G2 (SS) は, 検定統計量 G2 (CSS) と G2 (GS) の和に同等で
ある.
第 2 部 正方分割表における併合した表を用いた対称性に関する尺度
1. はじめに
行および列に対するカットポイント s および t (1 ! s < t ! r − 1) を用いて, カテ
$ %
ゴリを併合して r ×r 分割表を 3×3 分割表にする. 全部で r−1
(= (r −1)(r −2)/2)
2
枚の併合した分割表を作成することができ, それぞれの併合した分割表を Tst とす
(s,t)
る. このとき, Tst 表における (i, j) セル確率を Gij
(i = 1, 2, 3; j = 1, 2, 3) とする.
対称モデルが成り立たないとき, その隔たりの程度を測ることに関心があり, 第 2
部では併合した分割表を用いて対称性からの隔たりを測る尺度を提案した.
2. 提案する尺度
& & (s,t)
{pij + pji $= 0} を仮定する. ここに δst =
Gij (1 ! s < t ! r − 1) として,
i#=j
(s,t)
∗(s,t)
Gij
(s,t)
Gij
Gij
c(s,t)
=
, Gij = (s,t)
(i = 1, 2, 3; j = 1, 2, 3),
(s,t)
δst
Gij + Gji
とする. このとき, 対称モデルからの隔たりを測る尺度を次のように提案した:
1 ! ! (λ)
Ψ(λ) = $ r−1 %
Ψst (λ > −1),
2
ただし
(λ)
Ψst
1!s<t!r−1
")
(
* + 1 1 ,#
λ(λ + 1) ! ! ' ∗(s,t)
∗(s,t)
(λ)
c(s,t)
c(s,t)
= λ
Gij + Gji
Iij
Gij , Gji
;
,
,
2 −1
2 2
i<j
ここに
(λ)
Iij (·; ·) =
2
3
Gc(s,t) λ
1
Gc(s,t)
ij
λ(λ + 1)
ij
1/2
−1
c(s,t)
+ Gji
2
3
Gc(s,t) λ
ji
1/2
− 1 .
λ = 0 のときは λ → 0 の極限で定義した.
(λ)
部分尺度 Ψst (1 ! s < t ! r − 1) は Tst 表における対称性からの隔たりを測る尺
(λ)
c(s,t)
度である. ここに Iij (·; ·) は，{Gij
c(s,t)
, Gji } と {1/2, 1/2} の間の power-divergence
であり, 特に λ = 0 のときは Kullback-Leibler 情報量である. 提案尺度は以下の性
質をもつ.
• 0 ≤ Ψ(λ) ≤ 1
• Ψ(λ) = 0 ⇔ 対称モデルが成り立つ
• Ψ(λ) = 1 ⇔ 対称モデルからの隔たりが最大である
1.
Bowker, 1948
Stuart, 1955
McCullagh, 1978
Goodman, 1979
Agresti, 1983
1987
2
Tomizawa,
Tomizawa, 1993
Tomizawa, 1993
Miyamoto, Ohtsuka and Tomizawa, 2004
2
Tomizawa, Miyamoto, Yamamoto and Sugiyama, 2007
Miyamoto et al., 2004
Tomizawa et al., 2007
2
2.
R×R
(i, j)
pij
i = 1, . . . , R; j = 1, . . . , R
{Gij } i "= j
Gij =
i !
R
!
pst
(i < j),
Gij =
s=1 t=j
j
R !
!
s=i t=1
pst
(i > j).
2
-1
G+
i = ∆1
R
!
Ωk−i Gki
(i = 1, . . . , R − 1),
k=i+1
G+
i
R
!
=
Gik ,
k=i+1
2
-2
G−
i = ∆2
i−1
!
Λi−k Gki
(i = 2, . . . , R),
k=1
G−
i
=
i−1
!
Gik ,
k=1
∆ 1 = 1 ∆2 = 1
1
-1
1
-2
3.
2
1
t = 1, 2
2
2
-t
.
1
2
2
2
.
2
D-t
(
(
(
)
)
)
1.
McCullagh [3]
X1
(i, j)
(S)
pij = ψij
(PS)
X2
r r
pij (i = 1, . . . , r; j = 1, . . . , r)
(Bowker [2]; Bishop, Fienberg and Holland [1])
(i = 1, . . . , r; j = 1, . . . , r)
Gij =
i !
r
!
pst = P(X1 ≤ i, X2 ≥ j) (i < j),
r !
i
!
pst = P(X1 ≥ j, X2 ≤ i) (i < j).
s=1 t=j
Gji =
ψij = ψji
s=j t=1
S
Gij = Ψij
(i #= j)
m
Tahata and Tomizawa, [4])
αi
∆(m)
Ψij
i
αj−1
Gij =
Ψ
ij
pii = Ψii
Ψij = Ψji
(PS(m))
m (m = 1, . . . , r − 1)
(i < j),
pii = Ψii
Ψij = Ψji
(i > j),
(m)
∆i
=
m−1
%
k
∆ik .
k=0
PS(1)
(GPS)
PS(r − 1)
(McCullagh [3])
PS
2.
PS(m)
S
k
(MO(k))
k (≥ 1)
E(X1k ) = E(X2k ).
(Saigusa,
k=1
ME (
(CSS)
)
.
(Tomizawa, Miyamoto and Ouchi
[5])
(i = 1, . . . , r − 2).
Gi,i+2 = Gi+2,i
l
(MA(l))
l
(≥ 2)
1
2
µX
= µX
l
l ,
l = 2
VE (
l
t
µX
l = E((Xt − E(Xt )) ) (t = 1, 2).
)
.
(SE)
1
2
µX
µX
3
3
=
.
1 3/2
2 3/2
(µX
(µX
2 )
2 )
(KE)
1
2
µX
µX
4
4
= X2 2 .
1 2
(µX
(µ2 )
2 )
1
PS(m)
m (m = 1, . . . , r − 1)
,S
MO(k)
(k = 1, . . . , m)
1 (i) (ii) (iii)
.
(i) S
.
(ii)
m (m = 2, . . . , r − 1)
(l = 2, . . . , m)
CSS
(iii) r ≥ 5
PS(4)
ME
CSS
PS(m)
VE
CSS
ME
SE
MA(l)
KE
[1] Bishop, Y. M. M., Fienberg, S. E. and Holland, P. W. (1975). Discrete multivariate
analysis: theory and practice. MIT Press, Cambridge.
[2] Bowker, A. H. (1948). A test for symmetry in contingency tables. Journal of the American Statistical Association 43, 572-574.
[3] McCullagh, P. (1978). A class of parametric models for the analysis of square contingency
tables with ordered categories. Biometrika 65, 413-418.
[4] Saigusa, Y., Tahata, K. and Tomizawa, S. (2014). An extended asymmetry model for
square contingency tables with ordered categories. Model Assisted Statistics and Applications, to appear.
[5] Tomizawa, S., Miyamoto, N. and Ouchi, M. (2006). Decompositions of symmetry model
into marginal homogeneity and distance subsymmetry in square contingency tables with
ordered categories. Revstat: Statistical Journal 4, 153-161.