窓打切りされた再発事象のモデルとパラメータ の推定方法

中央大学大学院理工学研究科経営システム工学専攻 博士論文

窓打切りされた再発事象のモデルとパラメータ の推定方法

The models of recurrent events under the window censoring and their parameter estimation

阿部興

Ko ABE

学籍番号

14S7100002K

指導教員 鎌倉稔成 教授

1

目次

第

1

7

1.1

. . . . 7

1.2

. . . . 8

第

2

10 2.1

. . . . 10

2.2

. . . . 10

2.3

. . . . 11

2.4

. . . . 13

2.5

. . . . 16

2.6

. . . . 21

2.7

. . . . 24

2.8

. . . . 25

第

3

26 3.1

. . . . 26

3.2

. . . . 26

3.3

. . . . 29

3.4

. . . . 29

3.5

. . . . 31

3.6

. . . . 34

3.7

. . . . 42

第

4

50 4.1

. . . . 50

4.2

. . . . 50

4.3

. . . . 50

4.4

. . . . 51

4.5

. . . . 53

目次

2

4.6

. . . . 54

4.7

. . . . 54

4.8

. . . . 56

第

5

65 5.1

. . . . 65

5.2

. . . . 65

5.3

. . . . 67

5.4

. . . . 68

5.5

. . . . 71

5.6

. . . . 90

5.7

. . . . 91

第

6

93

95

97 A

Laslett (1982)

. . . . 97

B MTTR

. . . . 99

3

図目次

2.1

. . . . 12

2.2 2

. . . . 15

2.3

Z(t)

. . . . 18

3.1

I

. . . . 28

3.2

. . . . 29

3.3

. . . . 31

3.4

.

,

. . . . 33

3.5 Mauldon(1998)

. . . . 39

4.1

. X 1

X 2

.

. . . . . 51

4.2

.

. . . . 51

4.3

.

,

η 1 = 8, η 2 = 2

.

10

. . . . . 56

4.4

.

,

η 1 = 8, η ₂ = 2,m ₁ = 3,m ₂ = 3,

10

. . . . 57

4.5

.

,

η ₁ = 8, η ₂ = 2,m ₁ = 5,m ₂ = 5,

10

. . . . 57

4.6

.

,

η ₁ = 8, η 2 = 2,m 1 = 7,m 2 = 7,

10

. . . . 58

4.7

.

,

η 1 = 8, η 2 = 2,m 1 = 0.8,m 2 = 0.8,

10

. . . . 58

4.8

.

,

η 1 = 8, η 2 = 2,m 1 = 0.5, m 2 = 0.5,

10

. . . . 59

4.9

.

,

η 1 = 8, η 2 = 2,m 1 = 0.3,m 2 = 0.3,

10

. . . . 59

4.10

.

,

η ₁ = 8, η ₂ = 2,

5

. . . . . 60

図目次

4 4.11

.

,

η 1 = 8, η 2 = 2,

サンプルサイズ

100

. . . . 60

4.12

.

,

η 1 = 8, η 2 = 2,m 1 = 5,m 2 = 5,

5

. . . . 61

4.13

.

,

η 1 = 8, η 2 = 2,m 1 = 0.5,m 2 = 0.5,

5

. . . . 61

4.14

.

,

η 1 = 8, η ₂ = 2,m ₁ = 5,m ₂ = 5,

100

. . . . 62

4.15

.

,

η ₁ = 8, η ₂ = 2,m ₁ = 0.5,m ₂ = 0.5,

100

. . . . 62

4.16

.

,

k ₁ = 2, k 2 = 2,σ 1 = 4,σ 2 = 1,

5

. . . . 63

4.17

.

,

k 1 = 2, k 2 = 2,σ 1 = 4,σ 2 = 1,

10

. . . . 63

4.18

.

,

k 1 = 2, k 2 = 2,σ 1 = 4,σ 2 = 1,

100

. . . . 64

5.1 Rootz´ en & Zholud (2016)

5

. . . . . 66

5.2 µ 0

.

, Rootz´ en & Zhould

. . . . . 74

5.3 µ ₀

.

,

. . . . . 75

5.4

m

.

(m = 0.5)

, Rootz´ en & Zhould

. . . . . 75

5.5

η

.

m = 0.5

, Rootz´ en & Zhould

. . . . . 76

5.6

m

.

(m = 2)

, Rootz´ en & Zhould

. . . . . 76

5.7

η

.

m = 2

, Rootz´ en & Zhould

, . . . . 77

5.8

m

,

m = 0.5

,

. 78 5.9

η

.

m = 0.5

,

. 79 5.10

m

.

m = 2

,

. . . 79

5.11

η

.

m = 2

,

. . . 80

5.12

k

.

k = 3

, Rootz´ en & Zhould

. . . . 80

5.13

σ

.

k = 3

, Rootz´ en & Zhould

. . . . 81

図目次

5 5.14

k

.

k = 0.5

, Rootz´ en & Zhould

タイプの打切りの場合

. . . . 81 5.15

σ

.

k = 0.5

, Rootz´ en & Zhould

イプの打切りの場合

. . . . 82 5.16

k

.

k = 3

,

. . . . . 82 5.17

σ

.

k = 3

,

. . . . . 83 5.18

k

.

k = 0.5

,

場合

. . . . 83 5.19

σ

.

k = 0.5

,

場合

. . . . 84 5.20 MTTR

.

, Rootz´ en & Zhould

場合

. . . . 85 5.21 MTTR

.

(m = 2)

, Rootz´ en & Zhould

プの打切りの場合

. . . . 85 5.22 MTTR

.

(m = 0.5)

, Rootz´ en & Zhould

イプの打切りの場合

. . . . 86 5.23 MTTR

.

(k = 3)

, Rootz´ en & Zhould

の打切りの場合

. . . . 86 5.24 MTTR

.

(k = 0.5)

, Rootz´ en & Zhould

プの打切りの場合

, . . . . 87 5.25 MTTR

.

,

. . . . 87 5.26 MTTR

.

(m = 2)

,

. . . . . 88 5.27 MTTR

.

(m = 0.5)

,

. . . . 88 5.28 MTTR

.

(k = 3)

,

. . . . 89 5.29 MTTR

.

(k = 0.5)

,

. . . . 89 5.30

.

,

. . . . 92 B.1 MTTR

.

m = 0.5

, Rootz´ en & Zhould

イプの打切りの場合

. . . . . 100

B.2 MTTR

.

m = 0.5

,

. . . . 100

6

表目次

3.1

.

. . . . 41

3.2

.

. . . . 42

3.3

.

. . . . 43

3.4

RMSE

.

. . . . 44

3.5

RMSE

.

. . . . 45

3.6

RMSE

.

. . . . 46

3.7 µ

.

. . . . 47

3.8

.

. . . . 48

3.9

.

. . . . 49

5.1 Rootz´ en & Zholud

2016

. . . . . 78

5.2

AIC. . . . 91

A.1

(k = 2

) . . . . 97

A.2

(k = 0.5

) . . . . 98

A.3

(m = 2

) . . . . 98

A.4

(m = 0.5

) . . . . 99

7

第

1

序論

1.1

本研究の動機はヒートシールの強度の分析から得られた

.

,

.

,

,

,

.

,

,

.

,

.

,

.

.

,

Laslett’s line segment problem

Van Zwet, 2004

.

Laslett’s line segment problem

.

多く存在する

. Laslett (1982)

. Mauldon (1998)

, Svensson et al.(2006)

.

生存時間分析や信頼性工学などの分野では

,

,

,

.

,

,

.

,

,

,

.

,

.

.

.

(reccurent event)

.

観測が終了するまでにイベントの生起を確認できなかった場合

,

.

,

,

,

window-censored

(Zhao, 2006, Zhu et al., 2014,

Rootz´ en & Zholud, 2016).

,

第

1

8

1980

. Laslett (1982)

,

ン過程を仮定して

,

0

,

.

,

. Vardi (1982)

,

.

, Zhao (2006)

.

,

,

,

.

このように再生過程の研究は広くなされているが

,

.

2

.

,

.

.

,

,

2

(Barlow & Preschan,1965). Peuter (2013)

,

. Rootz´ en & Zholud (2016)

1,

0

,

0

.

1.2

第

2

,

.

,

,

,

.

,

.

第

3

.

Laslett(1982)

Mauldon(1998)

.

Laslett(1982)

,

.

第

4

(2016)

,

.

,

.

.

第

5

,

(2017)

, Rootz´ en & Zholud (2016)

,

(2016)

.

,

(2016)

Rootz´ en & Zholud (2016)

,

.

0

,

1

,

.

,

,

,

. Rootz´ en & Zholud (2016)

0

第

1

9

,

1

,

部・鎌倉

(2016)

,

1

,

0

,

1

布を同時に推定するアプローチを提案している

.

,

0

,

.

また窓打切り状況下での回帰型のモデルについては

, Zhu et al. (2014)

. Zhu et

al. (2014)

,

.

ルは阿部・鎌倉

(2016)

Rootz´ en & Zholud (2016)

,

.

(2017)

, 5.6

,

.

10

第

2

対象とする確率過程

2.1

この章では

,

,

2

,

,

,

.

,

.

2.2

,

,

.

ポアソン過程はイベントの生起の間隔が常に同一の指数分布に従う

.

,

,

.

.

.

,

, 2

,

,

.

この章の内容は

, Cox and Isham (1980), Ross (1995),

(1996), Beichelt & Paul (2002)

.

2.2

ここでは再発事象を考えるときに用いられる基本的な量である

,

.

確率密度関数

f (t)

F (t)

, 1 − F (t)

(survival function)

.

h(t)

,

h(t) = lim

∆t → 0

Pr(t ≤ T < t + ∆t | T ≥ t)

∆t

と定義される

.

t

,

第

2

11

. T

,

h(t) = lim

∆t → 0

1

∆t

Pr(t ≤ T < t + ∆t ∩ T ≥ t) Pr(T > t)

= lim

∆t → 0

Pr(t ≤ T < t + ∆)

∆t

/ Pr(T > t)

= f (t) 1 − F (t)

.

H (t) =

∫ t 0

h(u) du

.

,

1 − F(t) = exp( − H(t))

.

,

h(t) = f(t) S(t) = − d

dt log { 1 − F (t) }

,

,

exp (

−

∫ t 0

h(u) du )

= { 1 − F(t) }

.

以上の議論より

,

,

,

,

,

1

.

2.3

ポアソン過程は確率過程の中で最も基本的で重要なものである

.

N(t)

[0, t)

.

ても構わない

.

, N (t)

.

1.

,

.

s < t ≤ u < v

N (t) − N (s)

N (v) − N(u)

.

2.

,

.

P (N(t + h) − X (t)) = λh + o(h),

o(h)

h

, λ

. 3.

2

,

1

,

無視できるほど小さい

.

P(N (t + h) − N (t) ≥ 2) = o(h).

4.

0

,

0

.

N (0) = 0.

第

2

12

t

t

t

N(t)

図

2.1

上記の仮定にもとで

,

t ₁ , t ₂ , . . .

λ

.

1, 2, 3

,

P (N (t + h) = 0) = (1 − λh)P (N(t) = 0) + o(h)

.

.

P ^′ (N(t) = 0) = − λP (N (t) = 0)

,

P(N (t) = 0) = exp( − λt + C)

（

C

.

4

P(N (0) = 1)

, P (N (t) = 0) = exp( − λt)

.

. P (N (t) = 0)

P (t i > t)

, t 1

. 1

,

P (t 2 > t | t 1 > s) = P (N (t + s) − N(s) = 0 | N (s) = 1)

= P (N (t + s) − N(s) = 0)

= e ^{− λt}

.

.

図

2.1

t 1 , t 2 , . . . ,

N(t)

.

2.4

,

.

定理

1.

X

,

P(X > s + t | X > s) = P(X > t), s ≥ 0, t ≥ 0. (2.1)

,

.

まず

P (X > s + t | X > s) = P(X > s + t)/P (X > s) = P(X > t)

.

P (X > s + t | X > s) = P (X > t)

,

P (X > s + t | X > s) = P (X > s)P (X > t).

第

2

13 s = t = 0

,

P (X > 0) = 1.

, G(x) = P (X > x)

,

G(s + t) = G(s)G(t).

これを満たすためには

, G(t) = e ^{− λt}

.

次に

,

X 1 , X 2 , . . . ,

λ

N (t)

λt

. S _n = X ₁ + · · · + X _n , S ₀ = 0

. N (t) ≥ n

S _n ≤ t

.

N(t) = 0

. N (t) = 0

,

, t

,

Pr(N(t) = 0) = Pr(X 1 > t) = e ^λt

N (t) = n(n = 1, 2, . . .)

.

P (N(t) = n) = P (N (t) ≥ n) − P (N(t) ≥ n + 1)

= P (S n ≤ t) − P (S n+1 ≤ t) S n = X 1 + · · · + X n

λ, n

,

P (N(t) = n) =

∫ t 0

λ ⁿ x ^{n − 1}

(n − 1)! e ^{− λx} dx −

∫ t 0

λ ⁿ⁺¹ x ⁿ

n! e ^{− λx} dx

=

∫ t 0

λ ⁿ x ^{n − 1}

(n − 1)! e ^{− λx} − λ ⁿ⁺¹ x ⁿ

n! e ^{− λx} dx

= [ (λx) ⁿ

n! e ^{− λx} ] t

0

= (λt) ⁿ n! e ^{− λt}

λt

.

2.4

ポアソン過程ではイベントの生起の間隔は常に同一の指数分布に従う

.

,

,

.

非負の整数の値を取る連続時間の確率過程

{ X (t), t ≥ 0 }

.

s, t ≥ 0

i, j, x(u), 0 ≤ u ≤ s

,

Pr { X (t + s) = j | X (s) = i, X(u) = x(u), 0 ≤ u < s } = Pr { X (t + s) = j | X (s) = i } .

,

(continuous-time Markov chain)

.

,

,

. s

, t + s

,

.

加えて

, P { X(t + s) = j | X (s) = i }

s

,

(stationary)

(homogeneous)

.

第

2

14

.

,

[s, s + t]

,

t

留まる確率と一致する

.

i

τ i

,

s, t ≥ 0

.

Pr { τ i > s + t | τ i > s } = Pr(τ i > t).

従い

τ i

,

.

以下の式で表される

t

i

j

.

P ij (t) = P { X(t + s) = j | X(s) = i } . (2.2)

P ij (t)

,

,

.

v _i = lim

∆t → 0

1 − P _ii (∆t)

∆t .

q ij = lim

∆t → 0

P ij (∆t)

∆t , i ̸ = j.

P ij (t + h) =

∑ ∞ k=0

P ik (t)P kj (h)

,

P _ij (t + h) − P _ij (t) = ∑

k

P _ik (t)P _kj (h) − P _ij (t)

= ∑

k ̸ =j

P _ik (t)P _kj (h) − { 1 − P _jj (h) } P _ij (t)

,

lim

h → 0

P ij (t + h) − P ij (t)

h = lim

h → 0

 



∑

k ̸ =j

P ik (t) P kj (h)

h − 1 − P jj (h) h P ij (t)

 



和と極限の順序交換が可能であると仮定して

, P _ij ^′ (t) = ∑

k ̸ =j

q kj P ik (t) − v j P ij (t) (2.3)

(2.3)

(Kolmogorov’s forward equations)

.

コルモゴロフの前向き方程式を解くことにより

, P ij (t)

.

2

,P ij (t)

.

0

1

2

,

.

P ₀₀ ^′ (t) = µP 01 (t) − λP 00 (t)

= − (λ + µ)P 00 (t) + µ,

第

2

15

X(t)

1

0

t

2.2 2

ここで

2

P 01 (t) = 1 − P 00 (t)

.

1

,

e ^(λ+µ)t

,

.

e ^(λ+µ)t [P ₀₀ ^′ (t) − (λ + µ)P ₀₀ (t)] = µe ^(λ+µ)t

,

d

dt [e ^(λ+µ)t P ₀₀ ^′ (t)] = µe ^(λ+µ)t

,

.

e ^(λ+µ)t P ₀₀ ^′ (t) = µ

λ + µ e ^(λ+µ)t + C,

ここで

C

.

0

,P 00 (t) = 1

,C = λ/(λ + µ)

.

,

P 00 (t)

.

P ₀₀ ^′ (t) = µ

λ + µ + λ

λ + µ e ^{− (λ+µ)t} .

P ₁₁ (t)

.

P ₁₁ ^′ (t) = λ

λ + µ + µ

λ + µ e ^{− (λ+µ)t} .

2

.

0

,

1

.

2

2.2

.

十分に時間が経過した後

,

.

t lim →∞ P ₁₁ ^′ (t) = λ

λ + µ . (2.4)

この量は

,

,

.

4

.

第

2

16

2.5

2.3

,

,

.

,

,

,

.

.

{ X i , i = 1, 2, . . . , }

F(x)

. S 0 = 0, S n = ∑ n

i=1 X i

,

N(t) = sup { n | S _n ≤ t } (2.5)

とする

. N (t)

t

.

{ N (t), t ≥ 0 }

（

renewal process

. N (t)

M (t) = E[N (t)]

renewal function

. S n

F(t)

, S n

, F(x)

n

F ⁽ⁿ⁾ (t)

.

M (t)

,

M (t) =

∑ ∞ n=1

F ⁽ⁿ⁾ (t)

である

.

, I _{{ n,t }}

n

[0, t]

1,

0

,

.

E[N (t)] = E [ _∞

∑

n=1

I _{{ n,t }} ]

=

∑ ∞ n=1

E[I _{{ n,t }} ]

=

∑ ∞ n=1

Pr { I _{{ n,t }} = 1 }

=

∑ ∞ n=1

Pr { S n ≤ t }

=

∑ ∞ n=1

F ⁽ⁿ⁾ (t).

第

2

17

M (t)

,

M (t) =

∑ ∞ n=1

F ⁽ⁿ⁾ (t)

= F(t) +

∑ ∞ n=1

F ⁽ⁿ⁺¹⁾ (t)

= F(t) +

∑ ∞ n=1

F ∗ F ⁽ⁿ⁾ (t)

= F(t) + F ∗ [ _∞

∑

n=1

F (n)(t) ]

= F(t) + F ∗ M (t)

とも書ける

.

F ∗ G(t)

F(t)

G(t)

.

. M (t) = F (t) +

∫ t 0

M (t − x)dF (x) (2.6)

式

(2.6)

renewal equation

.

,

ラプラス変換によって簡潔になる

.

F (t)

F ^∗ (s) =

∫ _∞

0

e ^{− st} dF (t) (2.7)

と表すことにする

. M (t)

, M ^∗ (s) =

∑ ∞ n=1

[F ^∗ (s)] ⁿ

となる

.

,

. M ^∗ (s) = F ^∗ (s)

1 − F ^∗ (s)

F ^∗ (s)

,

F ^∗ (s) = 1 − µs + o(s)

,

s lim → 0 sM ^∗ (s) = F ^∗ (s)

{ 1 − F ^∗ (s) } /s = 1

µ (2.8)

となる

.

ラプラス変換について

,

(Feller, 1971; XIII

5

).

定理

2. F (t)

[0, ∞ )

. F (t)

F ^∗ (t)

.

α

,

s → lim +0 s ^α F ^∗ (s) = C (2.9)

第

2

18

Z(t)

0 S

t S

図

2.3

Z(t)

ならば

,

t lim →∞

F(t)

t ^α = C

Γ(α + 1) (2.10)

となる

.

定理

2

(2.8)

,

t lim →∞

M (t) t = 1

µ (2.11)

となる

. (2.11)

(elementary renewal theorem)

.

さらに

,

X _k

, (2.11)

,

(Ross, 1995; p.110).

t lim →∞ M (t + τ) − M (t) = τ

µ (2.12)

(2.12)

(Blackwell’s theorem)

.

,

t lim →∞

M (t + τ) − M (t)

τ = 1

µ

.

,

.

t lim →∞

dM (t) dt = 1

µ (2.13)

m(t) = dM (t)/dt

(renewal density)

,

t

.

2.5.1

打切りデータを扱うにあたって

,

(residual life)

.

,

Z(t) = S _{N (t)+1} − t (2.14)

とおくと

, Z (t)

t

(

2.3).

第

2

19 Z t

,

.

Pr(Z(t) > x) = 1 − F (t + x) +

∫ t 0

m(t − u)[1 − F (u + x)] du (2.15)

1

t

,

2

t

t − u

,

u + x

.

u ∈ [0, t)

.

t

,

.

t lim →∞ Pr(Z(t) > x) =

∫ _∞

x [1 − F(u)] du

µ . (2.16)

(2.16)

Ross(1996)

,

(equilibrium distribution)

.

2.5.2

指数分布の均衡分布

指数分布の均衡分布は簡単に求まる

.

µ

. F(x) = 1 − exp

(

− x µ

)

. (2.17)

均衡分布の生存関数は

,

G(x) =

∫ _∞

x

1 µ exp

(

− u µ )

du (2.18)

= exp (

− x µ

)

(2.19)

,

.

さらに逆も成り立つ

.

,

,

.

, F (x)

(0, ∞ )

F

µ

G(x) =

∫ _∞

x

1 − F (u)

µ du

とする．このとき，

1 − F (x) = G(x)

F(x) = 1 − e ⁻

(

)

. G(x)

d

dx G(x) = − 1 − F (x)

µ .

1 − F(x) = G(x)

d

dx G(x) = − 1 µ G(x).

これは

1

. G(x) = e ⁻

^x + C,

ここで

C

G(0) = 1 − F (0) = 1

C = 0

F (x) = 1 − e ⁻

1 − F (x) = G(x)

第

2

20

ワイブル分布の均衡分布

ワイブル分布の均衡分布の生存関数を求める

.

. F(x) = 1 − exp

{

− ( x

η ) m }

. (2.20)

ワイブル分布の平均は

ηΓ(1 + 1/m)

,

. G(x) =

∫ _∞

x

g(u) du = 1 ηΓ(1 + 1/m)

∫ _∞

x

exp {

− ( u

η ) m }

du. (2.21)

ここで

,

v = ( u

η ) m

とおくと

,

∫ _∞

x

exp {

− ( u

η ) m }

du = 1

m ηΓ(1/m, ( x η ) ^m )

,

.

G(x) = Γ

( 1 m ,

( x η

) m ) Γ ( ₁

m

) . (2.22)

ガンマ分布の均衡分布

次にガンマ分布の均衡分布の生存関数を求める

.

. F(x) = γ(k, x/σ)

Γ(k) (2.23)

ガンマ分布の平均は

kσ

,

g(x) = Γ(k, x/σ)

Γ(k + 1)σ (2.24)

ラプラス変換を用いると計算が容易になり

,

. G(x) =

∫ _∞

x

g(u) du =

∫ _∞

x

Γ(k, x/σ) Γ(k)kσ du

= 1

Γ(k + 1) (

Γ(k + 1, x/σ) − x

θ Γ(k, x/σ) )

(2.25)