POZNA MKA K VY POC TU BAYESOVSKE HO RIZIKA Ales LINKA TU Liberec, KPDM Abstrakt. V teto praci porovname dva bayesovske odhady funkce spolehlivosti v exponencialnm rozdelen z pohledu bayesovskeho rizika vypocteneho vzhledem k apriornmu gamma rozdelen ;(a; p) v situaci, kdy jsme v prpade jednoho z odhadu nepresne urcili hodnoty parametru apriornho rozdelen. Jako nastroj pro porovnan odhadu pouzijeme asymptotickou de ciency stanovenou na zaklade asymptotickych rozvoju pro bayesovske riziko uvazovanych odhadu. Abstract: In this paper we are interested in two Bayes estimators of reliability function in exponential distribution which have dierent a priori parameters. The asymptotic expansions of Bayes risk, computed with respect to the a priori distribution of gamma-type, are derived. For detailed comparison we use limit risk de ciency according to Lehmann (1983). Rezme: V to stat~e my zanimaems sravneniem dvuh baesovskih ocenok funkcii nadnosti v ksponencial~nom raspredelenii po vidu baesovskogo riska isqislenogo vzgldom k apriornomu gamma-raspredeleni. Baesovskovie ocenki otliqats vyborom apriornyh parametrov. Dl sravneni my ispol~zuem asimptotiqesku deficienci po Lemanu (1983).
1. Uvod
Uvazujme klasickou situaci, kdy vysledkem experimentu je uplny nahodny vyber X = (X1 ; : : : ; Xn )0 z exponencialnho rozdelen s hustotou 8 <
f (x; ) = :
1 exp ; x ; x > 0, 0; jinak,
(1)
kde 2 (0; 1) je neznamy parametr. Necht' c je kladne cslo. Jestlize X je nahodna velicina s hustotou (1), pak odpovdajc funkce spolehlivosti ma tvar
R(c) = R(c; ) = P (X > c) = exp f;c=g :
(2)
Dale predpokladejme, ze ztratova funkce ma tvar L(R(c); Rb) = [R(c) ; Rb]2 , kde R(c) oznacuje neznamou spolehlivost a Rb jej odhad. Riziko prslusne odhadu Rb, je-li skutecna hodnota spolehlivosti R(c), de nujeme vztahem
r R(c); Rb = E L R(c); Rb :
(3)
V tomto clanku budeme uvazovat bayesovsky odhad spolehlivosti a studovat jeho vlastnosti. Predpokladejme proto, ze parametr = ;1 je nahodna velicina s apriorn hustotou 8 <
q(; p; a) = :
ap p;1 e;a ; > 0; a > 0; p > 0; ;(p) 0; jinak.
(4)
Bayesovsky odhad R(c) dostaneme jako stredn hodnotu e; c vzhledem k aposteriornmu rozdelen, coz vede na odhad n+p nT + a n R R c) = nT + a + c : n b30
b30 (
(5)
Poznamenejme, ze oznacen prvnho indexu bylo zvoleno v souladu s oznacenm bayesovskeho odhadu v praci Hurt (1976) a v praci Antoch, Brzezina a Linka (1996), druhym indexem budeme rozlisovat ruzne bayesovske odhady. Bayesovske riziko odhadu Rb je de novano jako stredn hodnota rizika r vzhledem k apriornmu rozdelen q, v nasem prpade de novanem hustotou (4), tj.
% q; Rb = Eq r R(c); Rb =
Z Z
L R(c); Rb f (x; ) dx q(; p; a) d; (6) n
IR
kde f (x; ) je sdruzena hustota nahodneho vektoru (X1 ; : : : ; Xn )0 . 2. Asymptoticky rozvoj bayesovskeho rizika
V tomto odstavci budeme studovat chovan bayesovskeho odhadu vzhledem k bayesovskemu riziku. Pro dva ruzne bayesovske odhady urcme asymptoticke rozvoje pro jejich bayesovska rizika vzhledem ke kvadraticke ztratove funkci L a apriorn hustote (4). Oba odhady potom porovname pomoc asymptoticke de cience. Uvazujme nejprve situaci, kdy hledame asymptoticky rozvoj pro bayesovske riziko bayesovskeho odhadu Rb30 , pricemz toto bayesovske riziko poctame
vzhledem k apriornmu rozdelen de novanemu hustotou (4). Je-li X1 ; : : : ; Xn nahodny vyber z exponencialnho rozdelen s hustotou (1) a = ;1 , potom nahodna velicina n X 1 Tn = n Xi i=1
ma hustotu
8 <
h(t; ) = :
n nn tn;1 e; n t ; t > 0; > 0; ;(n) 0; jinak.
(7)
Ze znamych vlastnost gamma rozdelen a zrejmych upravach pro bayesovske riziko odhadu Rb30 dostavame
%(q; Rb30 ) = =
Z
1Z 1
0
0
a
= a +2c
p
n t + a n+p ; expf;c g2 h(t; )q(; p; a) d dt = nt+a +c
p ;(n + p) Z 1 xp;1 ;1 + a x n+p a n ; np ;(n) ;(p) 2 (n+p) dx: ( a + c ) x 0 1+
(8)
n
Konstrukce asymptotickeho rozvoje pro bayesovske riziko odhadu Rb30 funkce spolehlivosti vyuzva Taylorovu a Lebesquovu vetu. Protoze integrandy vystupujc v integralnch vyjadrench pro bayesovska rizika, jako funkce v 1=n, nejsou de novany v nule, tzn. v bode, ve kterem provadme konstrukci odhadu, je nutne odvodit modi kaci Taylorovy vety pro funkce spojite dode novatelne v tomto bode. Po nalezen rozvoje integrandu aplikujeme Lebesquovu vetu. Jedna se o postup technicky narocny a je nezbytne dokazat radu dlcch tvrzen. Vysledek je uveden v nasledujc vete.
Veta 1: Necht' apriorn hustota q je de novana vztahem (4), necht' n 2 IN a p; a; c 2 (0; 1). Polozme p2 %1 (q; Rb30 ) = a(ac+(p2 c+)p1)+2p ; %2 (q; R
b30 )
;
(9)
ap c2 p (1 + p) ;2a2 + 2c2 ; 2a2p ; 5c2p ; 4acp ; c2 p2 : (10) = 2(a + 2c)p+4
Potom pro asymptoticky rozvoj bayesovskeho rizika odhadu Rb30 plat b b %(q; Rb30 ) = %1 (q;nR30 ) + %2 (qn; 2R30 ) + Op;a;c (n;3 ); pro n ! 1:
Dukaz. Podrobny dukaz viz Antoch, Brzezina a Linka (1996).
(11)
2
Predpokladejme nyn, ze parametr je nahodna velicina s apriorn hustotou (4) s parametry p1 a a1 . Potom prslusny bayesovsky odhad funkce spolehlivosti R(c) ma tvar n+p1 nT + a n 1 R R c) = nT + a + c : n 1 b31
b31 (
(12)
To odpovda situaci, kdy jsme nepresne urcili apriorn parametry. Jako apriorn parametry parametry jsme zvolili p1 ; a1 , ale spravne jsme meli pouzt p; a. Nyn vypocteme asymptoticky rozvoj bayesovskeho rizika odhadu Rb31 vzhledem ke kvadraticke ztratove funkci L a apriorn hustote (4). Po kratkem vypoctu pro riziko odhadu Rb31 dostavame
%(q; Rb31 ) = =
Z
1Z 1
0
0
a
= a +2c 2
p
;
n t + a1 n+p1 ; expf;c g2 h(t; )q(; p; a) d dt = n t + a1 + c
p ;(n + p) Z 1 xp;1 ;1 + a1 x n+p1 a n + np ;(n) ;(p) 2 (n+p1 ) ( a 0 1 +c) x 1+ x n+p1
a1 n (n+p) ( a 1 + nx
64 1 +
n
n+p1 3
;2
7
1 + (a1 +nc) x
(n+p) 5 dx: (( a + c ) x 1+ n
(13)
Pro odvozen asymptotickeho rozvoje pouzijeme podobneho postupu jako jsme pouzili pro odhad Rb30 . Dostaneme nasledujc tvrzen.
Veta 2 Necht' apriorn hustota q je de novana vztahem (4), necht' n 2 IN a p; a; c 2 (0; 1). Polozme p 2 p %1 (q; Rb31 ) = a(ac+(p2c+)p1) +2 ; %2 (q; Rb31 ) =
;
= 12a2 + 6a4 ; 12a1 2 ; 12a2a1 2 + 6a1 4 + 48ac + 48a3c ; 48a1 c; ; 48a2 a1 c ; 48aa1 2 c + 48a1 3 c + 88a2c2 ; 192aa1 c2 + 96a1 2 c2 + + 8c4 + 10a2p + 5a4 p ; 10a1 2 p ; 10a2a1 2 p + 5a1 4 p + 40acp + + 40a3 cp ; 40a1 cp ; 40a3cp ; 40a1 cp ; 40a2a1 cp ; 40aa1 2 cp + + 40a1 3 cp + 40a2c2 p ; 128aa1 c2 p + 80a1 2 c2 p ; 80ac3p + 64a1 c3 p ; ; 20c4 p + 2a2p2 + a4 p2 ; 2a1 2 p2 ; 2a2 a1 2 p2 + a1 4 p2 + 8acp2 + + 8a3 cp2 ; 8a1 cp2 ; 8a2a1 cp2 ; 8aa1 2 cp2 + 8a1 3 cp2 + 8a2c2 p2 ; ; 16aa1 c2 p2 + 16a1 2 c2 p2 + 32a1 c3 p2 + 28c4p2 + 32a2c2 p1 ; ; 32aa1 c2 p1 + 64ac3p1 ; 64a1 c3 p1 ; 16aa1 c2 pp1 ; 32ac3pp1 ; ; 32a1 c3 pp1 ; 64c4 pp1 + 8a2 c2 p1 2 + 8a2 c2 p1 2 + 32ac3p1 2 + + 32c4 p1 2 ap (a + 2c);4;p p (1 + p) : (14) Potom asymptoticky rozvoj bayesovskeho rizika odhadu Rb31 je dan vztahem b b %(q; Rb31 ) = %1 (q;nR31 ) + %2 (qn; 2R31 ) + Op;a;c (n;3 ); pro n ! 1:
Dukaz. Detailn odvozen viz Linka (1996).
(15)
2
3. Vypocet asymptoticke deficience
Vzhledem k tomu, ze odhady Rb30 a Rb31 jsou tzv. asymptoticky silne e cientn vzhledem k stredn kvadraticke odchylce, asymptoticke rozvoje bayesovskeho rizika pro odhady Rb30 a Rb31 maj tvar
%(q; Rb3i ) = an + bn32i + op;a;c(n;2 ); i = 0; 1; tj. koe cient u 1=n je pro oba odhady stejny. Pro detailn porovnan odhadu Rb30 a Rb31 muzeme uzt de cienci, blze
viz Lehmann (1983). Zhruba receno, de cience spoctena pro jistou pevnou dvojici odhadu bude v nasem prpade ukazovat o kolik vce (nebo mene) pozorovan vyzaduje odhad B , ma-li mt stejne bayesovske riziko jako odhad A zalozeny na vyberu rozsahu n. V praxi se obvykle uzva asymptoticka de cience pro n ! 1. Jestlize oznacme %n (q; A) a %n (q; B ) bayesovska rizika odhadu A a B , a plat-li %n (q; A) = nar + nrb+1 + o n;(r+1)
(16)
a
(17) %n (q; B ) = nar + nrc+1 + o n;(r+1) ; pak asymptoticka de cience odhadu B vzhledem k odhadu A je de novana
vztahem
Veta 3 Necht' p; a; c 2 (0; 1).
c;b dBA % = ar :
(18)
Pro asymptotickou de cienci odhad u Rb30 a R vzhledem k bayesovskemu riziku plat b31
dR%b31 Rb30 (p; a; p1 ; a1 ; c) = ;
= ;12a2 ; 6a4 + 12a1 2 + 12a2a1 2 ; 6a1 4 ; 48ac ; 48a3c + 48a2a1 c + + 48a1 c + 48aa1 2 c ; 48a1 3 c ; 96a2 c2 + 192aa1 c2 ; 96a1 2 c2 ; 10a2p ; ; 5a4 p + 10a1 2 p + 10a2a1 2 p ; 5a1 4 p ; 40acp ; 40a3cp + 40a1 cp + + 40a2a1 cp + 40aa1 2 cp ; 40a1 3 cp ; 48a2c2 p + 128aa1 c2 p ; ; 80a1 2 c2 p + 64ac3p ; 64a1 c3 p ; 2a2 p2 ; a4 p2 + 2a1 2 p2 + + 2a2 a1 2 p2 ; a1 4 p2 ; 8acp2 ; 8a3 cp2 + 8a1 cp2 + 8a2a1 cp2 + + 8aa1 2 cp2 ; 8a1 3 cp2 ; 8a2c2 p2 + 16aa1 c2 p2 ; 16a1 2 c2 p2 ; ; 32a1 c3 p2 ; 32c4p2 ; 32a2c2 p1 + 32aa1 c2 p1 ; 64ac3p1 + + 64a1 c3 p1 + 16aa1 c2 pp1 + 32ac3pp1 + 32a1 c3 pp1 + 64c4pp1 ; 8a2c2 p1 2 ; 32ac3p1 2 ; 32c4p1 2 8;1 c;2 (a + 2 c);2 (19)
Dukaz.
Vysledek dostaneme dosazenm do vzorce (18) podle vet 1 a 2. 2 4. Prklad
Prklad1.
Bott a Hass (1978) uvad doby do poruchy vstupnch tesncch zaklopnych ventilu pro jaderne reaktory. Kombinac techto historickych dat a rostouc urovne zatzen byly stanoveny pozadovane hodnoty pro 5 % a 95 % kvantil apriornho rozdelen pro intenzitu poruch. Pro hodnotu 5 % kvantilu byla stanovena hodnota 1:4 10;5 (poruch za hodinu) a pro 95 % kvantil hodnota 4:9 10;5. V monogra i Martz a Waller (1980) v kapitole 6 muzeme nalezt metodu, kterou tito autori vypracovali, pro stanoven parametru apriornho gamma rozdelen. Na zaklade tohoto postupu stanovme apriorn rozdelen jako gamma rozdelen s hustotou (4) s parametry a = 285 714 a p = 8:5.
Vzhledem k vyse uvedenemu rozdelen budeme uvazovat odhad Rb30 , ktery nam predstavuje odhad se spravne zvolenymi apriornmi parametry. Nyn porovname odhad Rb30 s odhadem Rb31 , kdy jsme se netre li presne do apriornch parametru. K porovnan pouzijeme asymptotickou de cienci (19). Na obrazcch 1{5 jsou znazorneny grafy de cience dR%b31 Rb30 (p; a; p1 ; a1 ; c) pro a = 285 714; p = 8:5 a c = 10 000; 60 000; 110 000; 160 000; 210 000. Pro tyto hodnoty c jsou rovnez uvedeny hodnoty de cience v krajnch bodech intervalu (0:7 p; 1:3 p) (0:7 a; 1:3 a). 2
20 350000
10 0 6
300000 a1 8
250000
p1 10 200000
Obr.1 dR%b
R (8:5; 285 714; p1; a1 ; 210 000)
31 b30
Tab.1 dR%b
R (8:5; 285 714; p1; a1 ; 210 000) pro vybrane hodnoty p1 a a1 .
31 b30
a1 p1
dbR% 31 bR30
0.7 a 0.7 p
0.7 a 1.3 p
1.3 a 0.7 p
1.3 a 1.3 p
8.87 21.88 24.32 11.31
30 20
350000
10 0 6
300000 a1 8
250000
p1 10 200000
Obr.2 dR%b
R (8:5; 285 714; p1; a1 ; 160 000)
31 b30
Tab.2 dR%b
R (8:5; 285 714; p1; a1 ; 160 000) pro vybrane hodnoty p1 a a1
31 b30
a1
0.7 a 0.7 p
p1
0.7 a 1.3 p
1.3 a 0.7 p
1.3 a 1.3 p
13.89 29.04 33.95 18.79
b31 R b30 dR %
60 40 350000 20 0 6
300000 a1 8
250000
p1 10 200000
Obr.3 dR%b
R (8:5; 285 714; p1; a1 ; 110 000)
31 b30
Tab.3 dR%b
R (8:5; 285 714; p1; a1 ; 110 000) pro vybrane hodnoty p1 a a1 .
31 b30
a1 p1
dbR% 31 bR30
0.7 a 0.7 p
0.7 a 1.3 p
1.3 a 0.7 p
1.3 a 1.3 p
28.13 46.29 58.71 40.56
150 100
350000
50 0 6
300000 a1 8
250000
p1 10 200000
Obr.4 dR%b
R (8:5; 285 714; p1; a1 ; 60 000)
31 b30
Tab.4 dR%b
R (8:5; 285 714; p1; a1 ; 60 000) pro vybrane hodnoty p1 a a1 .
31 b30
0.7 a 0.7 p
a1 p1
dbR% 31 bR30
0.7 a 1.3 p
1.3 a 0.7 p
1.3 a 1.3 p
93.75 116.37 168.44 145.81
4000 350000 2000 0 6
300000 a1 8
250000
p1 10 200000
Obr.5 dR%b
R (8:5; 285 714; p1; a1 ; 10 000)
31 b30
Tab.5 dR%b
R (8:5; 285 714; p1; a1 ; 10 000) pro vybrane hodnoty p1 a a1
31 b30
a1 p1
0.7 a 0.7 p
dbR% 31 bR30
3270.71
0.7 a 1.3 p
1.3 a 0.7 p
1.3 a 1.3 p
3300.74 5788.65 5758.02
5. Zaver
Z uvedenych obrazku 1{5 vyplyva, ze hlediska asymptoticke de cience ma bayesovsky odhad Rb31 funkce spolehlivosti R(c), tj. odhad, kdy jsme netre li apriorn parametry, mnohem hors chovan s klesajc hodnotou doby do poruchy c. V prpade c = 10000 a prehodnotme-li parametry o 30 % tato de cience cn dokonce 5758. Naopak pro velke hodnoty c se rozdly odhadu Rb30 a Rb31 straj. Literatura
[1] Antoch J., Brzezina M., Linka A., Asymptotic approximation of Bayes risk of estimators of reliability for exponentially distributed data, Statistics & Decision, (1996), to appear. [2] Bott T. F., Haas P. M., Initial Data Collection Eorts of CREDO : Sodium Value Failers, NCSR R20, (1978), National Center of Systems Reliability. [3] Hurt J., On estimation of reliability in exponential case, Aplikace matematiky 21 (1976), 263-272. [4] Linka A. Notice on Bayes estimators, Technical Report N.7, Technical University of Liberec, (1996). [5] Lehmann E. L., Theory of Point Estimation, John Wiley & Sons (1983), New York. [6] Martz H. F., Waller R. A., Bayesian Reliability Analysis, John Wiley & Sons (1982), New York.