Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Estimation with multivariate extreme value distributions with applications to environmental data Chen, Youshi 1991-12-31

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
UBC_1991_A6_7 C42_5.pdf [ 3.38MB ]
Metadata
JSON: 1.0302381.json
JSON-LD: 1.0302381+ld.json
RDF/XML (Pretty): 1.0302381.xml
RDF/JSON: 1.0302381+rdf.json
Turtle: 1.0302381+rdf-turtle.txt
N-Triples: 1.0302381+rdf-ntriples.txt
Original Record: 1.0302381 +original-record.json
Full Text
1.0302381.txt
Citation
1.0302381.ris

Full Text

ESTIMATION W I T H MULTIVARIATE E X T R E M E V A L U E DISTRIBUTIONS, WITH APPLICATIONS TO ENVIRONMENTAL DATA By YOUSHI C H E N B.S c , Beijing Institute of Technology, 1983 M.S c , Beijing Institute of Technology, 1986  A THESIS S U B M I T E D IN P A R T I A L F U L F I L L M E N T O F THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE in T H E F A C U L T Y O F G R A D U A T E STUDIES D E P A R T M E N T O F STATISTICS  We accept this thesis as conforming to the required standard  T H E UNIVERSITY OF BRITISH COLUMBIA September 1991 ©Youshi Chen, 1991  In presenting  this thesis in partial fulfilment of the requirements for an advanced  degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may department  or by  his or her representatives.  be granted by the head of my  It is understood that  copying or  publication of this thesis for financial gain shall not be allowed without my permission.  Department of The University of British Columbia Vancouver, Canada Date  DE-6 (2/88)  written  Abstract Several parametric families of multivariate extreme value distributions (Hiisler and Reiss 1989, Tawn 1990, Joe 1990a, 1990b) have been proposed recently. Applications to multivariate extreme value data sets are needed to assess the adequacy of the known families in their fit to data. Different families are compared in their range of multivariate dependence and their ease of use for maximum likelihood estimation.  Some useful  conclusions have been made from experience with several environmental data sets.  n  TABLES OF CONTENTS Abstract  ii  Table of Contents  iii  List of Tables  iv  Acknowledgement  v  Chapter 1. Introduction  1  Chapter 2. Multivariate Extreme Value Distributions  3  2.1 Introduction  3  2.2 Multidimensional Distributions  3  2.3 Copulas  6  2.4 Univariate Extremes and the Generalized Extreme Value Distribution  7  2.5 Multvariate Extreme Value Distributions  11  2.6 Some Relation Between Extreme Value Theory and Central Limit Theory 18 2.7 Parametric Families of Multivariate Extreme Value Distributions  20  2.8 Estimation for Multivariate Extreme Value Distribution  23  Chapter 3. Maximum likelihood Estimation With Husler-Reiss Model 25 3.1 Bivariate Case  26  3.2 Trivariate Case  28  3.3 Four Dimensional Case  30  3.4 Computer Implementation of the Husler-Reiss model  32  3.5 Data Analysis  32  Chapter 4. Analysis with Missing Data  44  Chapter 5. Conclusion and Further Research  51  Bibliography  52  in  List of tables Table 1. Bay Area-Ozone. (n=120)  38  Table 2. GVRD-Sulphur Dioxide. (n=70)  39  Table 3. GVRD-Nitrogen Dioxide. (n=84)  40  Table 4. GVRD-Ozone. (n=84)  41  Table 5. GVRD-Station 5 (Confederation Park). (n=64)  42  Table 6. GVRD-Station 7 (Anmore). (n=76)  42  Table 7. GVRD-Station 9 (Rocky Point Park). (n=63)  43  Table 8. GVRD-Sulphur Dioxide. (n=104)  48  Table 9. GVRD-Nitrogen Dioxide. (n=84)  49  Table 10. GVRD-Ozone. (n=84)  ...50  iv  ACKNOWLEDGMENT I would like to express my deep gratitude to my supervisor, Professor Harry Joe, for his constant encouragement, guidance and assistance in my study and in the producing of this thesis. Very special thanks are due Professor Mohan Delampady for reading the manuscripts and making constructive comments. Finally, I would like to thank the Department of Statistics, The University of British Columbia, for supporting me with financial aid and encouragement.  v  Chapter 1 Introduction Extreme value theory has a great number of applications. We list two cases where the largest or the smallest "measurement" are of interest. 1. A i r Pollution. A i r pollution concentration is expressed in terms of proportion of a specific pollutant in the air. Concentrations are recorded at equal time intervals, and it is required by law to keep the largest annual concentrations below given limits. 2. Natural Disasters. Floods, heavy rains, extreme temperatures, extreme atmospheric pressures, winds and other phenomena can cause extensive human and material loss. Communities can take preventive action to minimize their effects even if such disasters cannot be completely avoided. In dams, dikes, canals, and other structures the choice of building materials and methods of architecture can take some of these potential disasters into account. Engineering decisions that confront such problems should be based on a very accurate theory, because inaccuracies can be very expensive. For example, a dam built at a huge expense may not last long before collapsing. In example 1, at each time point, there may be a vector X of measurements, consisting of concentrations of a pollutant at several air quality monitoring stations in a network or consisting of concentrations of several pollutants at a given monitoring stations. In example 2, there may be a vector X of measurements at several positions in a spatial grid. In the examples, we have over time a number n of random measurements X i , ' . . . , X „ , and the behaviour of either Z = m a x ( X i , . . . , X ) or W n  n  n  = m i n ( X i , . . . , X ) can be of n  interest. Here max and min refer to componentwise maxima and minima; see Chapter 2 for details. Statistical inferences of the upper or lower tail of the multivariate distribution of X can be made from using extreme value theory to model the multivariate distribution of  1  Z or W . n  n  There has been much recent research in multivariate extreme value theory; see Galambos (1987, Chapter 5) and Resnick (1987, Chapter 5) for some general theory. The set of possible of multivariate extreme value distributions has been shown to be of infinite dimension. Smith, Tawn and Yuen (1990) study a density estimation based method for nonparametrically fitting a multivariate extreme value distribution. Except in the bivariate case, the method was not too promising. A n alternative is to find a finite-dimensional parametric subfamily of the set of multivariate extreme value distributions. A good subfamily would be "dense" in the complete family. Parametric families have recently been proposed by Hiisler and Reiss (1989), Tawn (1990), Joe (1990a,b), Coles and Tawn (1991). Experience is now needed to assess the adequacy of fit of these families to multivariate extreme value data sets. Joe (1990b) did some work on this but did not include the Hiisler and Reiss family in his comparisons. One goal of this thesis is to do maximum likelihood estimation with the Hiisler and Reiss family. This turns out be very difficult and details are given in Chapter 3, where comparisons are also made with results in Joe (1990b). Another goal of the thesis is to deal with missing values in computing maxima. Joe (1990b) deleted cases where there were many missing values so that maxima could not be properly computed. These maxima were assumed to be missing at random. Here, in Chapter 4, we treat these cases as right-censored maxima, that is the actual maxima is known to be above the maxima of all non-missing observations. The thesis proceeds next (in Chapter 2) with definitions and properties of general multivariate distributions and of multivariate extreme value distribution. The last chapter (Chapter 5) of this thesis makes some final conclusions and points out a direction for further research.  2  Chapter 2. Multivariate Extreme Value Distribution 2.1 Introduction In this chapter we describe some general properties of multidimensional distributions and copulas. Then we introduce some univariate extreme value theory and the generalized extreme value distribution. Following this, we concentrate on multivariate extreme value distributions, the main topic of this chapter. 2.2 Multidimensional Distributions Let us define the p-dimensional random variable X as the vector X = (Xi,... The distribution function -F(x) = F(xi,...,  P  x ) is defined as p  F(x) = P ( X < x) = P(X where X < Y means Xi < Yi,  ,X )  <x ...,X <  X  u  x ),  p  p  for 1 < i < p.  For future reference, we also define X.+ Y = (Xi + Y ,...,X  p  X Y = (X\Yi,...,  p  1  + Y ), p  X Yp),  and X/Y  =  (X /Y ,...,X /Y ). 1  1  P  P  Random variables with the same distribution as X will be denoted by X i , . . . , X , p  and the components of X j by Xij, that is, Xij is the ith component of X j , or X j = {Xij, • • • ,X j). p  The order statistics of the ith. component (Xn,... ,Xi ) n  x<£ <xW  <...<xj \ n  and we denote Wi — XJ;l\  Z{ = 3  X$.  are  In extreme value theory, an objective is to investigate the existence of the asymptotic distribution of W  =  n  For  (W ,... W ) l  1  Z  p  — (Zi, • • • > Zp).  n  this, we need some theory of multidimensional distribution functions.  A multidimensional distribution function F ( x ) of a random vector X has the following  elementary properties: PI  (Bounds): 0 < F ( x i , . . . , x ) < 1  Vxi,..., x.  P2  (Monotonicity): F is nondecreasing in each of its arguments x,-, 1 < i < p.  p  P3 (Limits): linx^-oo F ( x  p  ...,x ) = 0. p  u  P4 (Marginal distributions): If Xj —• +oo, then F ( x ) tends to an (p—l)-dimensional distribution, which is the distribution of the vector obtained by removing its jth component. We can obtain the univariate marginal distribution FJ(XJ) by letting each x^, i ^ j tend to +oo. P5 (Density): If F has derivatives of order p, then f  =  >0,  dx\ • • •v dx.  and / is the probability density function. P6 (Rectangles):  P(a <X <b ,...,a <X <b ) 1  l  1  p  where for a subset I of { 1 , . . .  p  p  = F(b ,b ,...,b )-J2Fi 1  2  ^F  i j  p  p  + ... +  i=l  + (-lYF(a ,...,a )>0, 1  p  F j is the value of F ( s i , . . . , Sk) with s; = a; if i 6 / ,  and Si = b{ if i £ I. P7 (the Frechet Bounds): p  max(0, J2 i( i) F  x  ~ P + !) < ( ^ F x  •••  x ) < min(Fi(a;i),..., p  F (x )). p  p  If Fi,...,F  are all continuous, then the upper bound corresponds to distribution of  p  (Xi, F ~ F i ( X i ) , . . . , F ^ F ^ X i ) ) . In general the lower bound is not a multivariate dis1  2  tribution function for p > 3. For example, suppose F is a three-dimensional distribution function with marginal functions satisfying Fi(xi)  F (x ) 3  = 1/2,  ^2(^2)  = 1/2 and  = 1/2. Then  3  F(+oo,  +00,  - Fi(ari) - F2(x2)  +00)  + max(0, F ^ ) max(0, F2{x2) =  + F {x ) 2  2  + F (x ) 3  1-3/2 =-1/2  3  -  F (x ) 2  3  - 1) + max(0, F i ( x i ) + F3{x3) - 1) - max(0, F f a ) + F2(x2)  - 1) +  + F (x ) 3  3  - 2)  <0,  which contradicts Property 6. P8  (Survival function): Let Bj be the event (Xj  > Xj) and j(k)  = (ji,...  ,jk) be a  subset of { 1 , . . . ,p}, then let  G  j ( k )  (x ,..., h  X j k  ) =P(B  n • • • n B ),  h  5 (x)  = 1,  0  jk  and  k(x)  S  For  =  G  i<ji<h-ik<p  h ( h • •• >  p > 2, then  F(x ,...,x ) 1  =  p  x  p  >  x  < <Pk  3k)>  1  J2(- ) k(xi,...,x ). l ks  p  k=0  In addition, for any integer 0 < q < (p — l ) / 2 ,  E ( - i ) ^ W < f(x) < E ( - i ) ^ ) . f c  f c  fc=0  k-0  The proof of most of the above results is straightforward; for the proof of properties P7  and P8, readers are referred to Galambos (1987).  5  2.3 Copulas Assume that F(x) is a p-dimensional distribution function with univariate marginals Fi(x{), for 1 < i< p. Let C(y) be a p-dimensional function over the unit cube 0 < yi < 1, 1 < ^ < P, and such that it increases in each of its variables and F(x ,...,x ) 1  =  p  C[F (x ),...,F (x )]. 1  1  p  p  Then the function C(y) is called a copula of JP(X). If C is a copula, then G(Vu•••,y ) P  =  C[G (y ),...,G (y )] 1  1  p  p  is a multivariate distribution with univariate margins G,-, for i — 1, ...,p, where Gi is an arbitrary univariate distribution function. When needed, we write CF = CF^Y) = C(y) for ^(x). Note that some people call C(y) a dependence function. If F(x\,...,  x ) is a continuous p-variate cumulative distribution function with unip  variate margins Fj(xj),j = 1,... ,p, then the associated copula is C( ,..., UL  tip) = F(Fr (u ),...,  F -\u )),  1  1  p  0 < m < 1, i = 1,... ,p.  p  This is a multivariate distribution with uniform (0,1) margins. Returning to extreme values, the ith. marginal of F (x) is F^Xi), n  F ( ,...,  xp) = CFn[F?( ,...,  n  Xl  F?(x )].  Xl  p  On the other hand, we have F (x ,..., p)  =  n  1  X  C [F(x ),...,F(xp)}. F  1  A comparison of these last two equations leads to CF4y)  =  c [yl F  6  / n  ,...,yl  / n  }.  and thus  2.4 Univariate Extremes and the Generalized Extreme Value Distribution Let  Xi,...,X  denote independent and identically distributed (iid) random vari-  n  ables, with cumulative distribution function F. Let = max(Xi,..., X ),  Z  n  By  W  n  n  = min(X ,..., X ). a  n  the iid assumption, we have  H (x)  = P(Z < x ) =  n  F {x) n  n  and  L (x)  = P(W  n  We  seek conditions on F(x)  {d } n  < x) = 1 - (1 -  n  F(x)) . n  to guarantee the existence of sequences { a } , {&„}, {c }, n  n  of constants such that, as n —> oo  lim H J a  n  n—*oo  + b x) =  H(x)  n  and  lim L (c n  exist for all continuity points of H ( x )  + d x)  n  = L(x)  n  and L(x)  respectively, where H ( x )  nondegenerate distribution functions. Equivalently, {a }, n  that (Zn  — a )/b n  n  and (Wn — cn)/dn  {&„},  and L(x)  {c } and { d } n  n  are  are such  converge in distribution as n —> oo.  A sufficient condition (see Galambos, 1987) is given next. Assume that there are sequences a , c , b > 0, and d > 0, of real numbers such n  n  n  n  that, for all x and y, n[l — F(an  lim  + b x)] = u(x), n  n—>-+oo  lim  nF(c  + d y)  n  n  = w(y)  n—*-+oo  exist. Then  lim  P(Z  n  < a + b x) n  n  = exp[—u(x)]  and  Jim P(W  - 1 -  < c + d y)  n  n  n  exp[-w(y)].  Example 1. ( The Exponential Distribution). Let Xi,... ,X  n  be iid random variables  with common distribution function  F(x)  = 1 - e~ ,  x > 0.  x  Then  1 - F{a  + b x) = e - e a n  n  h n X  n  .  In order to satisfy  n[l — F(a  lim  + b x)\ — u(x),  n  n  n—•+oo  we can choose a = logn. and b = 1. Hence, u(z) n  n  lim On  P(Z  and thus,  = e~ , z  < logre+ z) = exp(—e~ ). z  n  the other hand, we can choose c = 0 and d = 1/re, such that n  nF(c  lim  -f- c? y) =  n  n  n—>+oo  n  lim nFiyln)  lim n ( l — e~y^n)  =  n—t-+oo  = y.  n—t-+oo  Hence, lim  P ( W < c + d y) = n  n  n  n-+oo .  lim P { W n n—*+oo  < y/n) = 1 -  e~ . y  -  Example 2. Let X i , . . . , X be iid random variables with common distribution function n  F(x) To  = \-x~ , x > l . 1  determine the limiting distribution of Z , we note that with a = 0 and b = re, we n  n  obtain lim  nil - F(a  n  + b y)\ = lim re — n  = -,  y > 0.  n  and thus lim P ( Z  < ny) = e x p ( - - ) .  n  n-++oo  y  For maxima, it is well-known that there are only three types of nondegenerate distributions H(x) that can be univariate limiting extreme value distributions. These results are contained in books devoted to extreme value theory, such as Galambos (1987). The three types are: H (x) = <  x > 0,7 > 0,  exp(—x~)  lf1  0 H (x) = <  otherwise, x < 0,7 < 0,  exp(—(—x) i )  2n  x > 0,  1  H3,o(x) - exp(—e )  — oo < x <  x  + 0 0 , 7 = 0.  After adding location and scale parameters the above three types of distributions can be summarized together as the generalized extreme value distribution ( G E V ) : , H(x) = exp{-[max(l +  ^ ~  7  ^ , 0 ) ] ^ }  - 00 < 7 < + 0 0 .  ( 2 - 4 - 1 )  The parameter 7 can be interpreted as a tail index parameter of F. It measures the thickness of the tail of F; 7 is larger if F has a thicker tail, and vice versa. The tail thickness is related to the rate of convergence of 1 — F to 0 as x —> 00. For example, let F(x) =. [max(l + 7 s ,  1 -  0)]  _ 1 / 7  .  For 72 < 71 < 0 then — 1 / 7 1 > — 1 / 7 2 > 0 and  limi/T4  ———;— = OO. [max(l + 72X, O ) ] - ^  [max(l+7 x,0)]- ^ 1  1  1  2  For 72 < 0, 71 = 0, then li  m  ^-1/72  7  „, ,,— = 0 0 .  T,  [max(l +  x  7 x,0)]- /i'2 1  2  9  For 72 = 0, 71 > 0, then [max(l-f :r,0)]~  r  x  lim ^oo  e  e  1/71  7l  r  e  x  x  - hm . x^oo [max(l + 7 i x , 0 ) ] / '  -x  1  > lim  r—,— = oo. (x7i) /'  71  1  yi  For 71 > 72 > 0, then 0 > - I / 7 1 > - I / 7 2 , I/72 — I/71 > 0 and ,. h  [max(l + 7 i ^ , 0 ) ] " (1 + 7 i * ) ~ T\ ,,— = l i m — —-;— [max(l + 72Z, O)]- /^ (1 + 7 z)- /'y2 1 / 7 1  m  1 / 7 1  r  7  1  1  2  > lim i  1 +  l x X ) X  ' T = lim (1 + ; r )  1 / 7 2  7 l  a;-+oo M +7lX) l^  1  X—KX>  ~  = 00.  1 / 7 1  '  K  So, the tail is heavier as 7 increases. For 7 > 0, from Theorem 2.1.1 of Galambos (1987), we find a  n  - [ ( l / r a ) " - 1] so that 1  7  = 0 and b =  n  7  = 1 - F(b x)  F(b x) n  = [1 + ( n - l)x}7  n  ~ n- x- .  ih  x  l h  and  F { b x ) = (1 - F(b x)) n  n  = (1 - n " ^ - ^ ) " ™  exp(-a;- ).  1  n  n  1/7  For 7 < 0, from Theorem 2.1.2 of Galambos (1987), we find a = — 1 / 7 , n  —7 n _ 1  7  b = n  so that  F(x)  = (1 + 7 * ) "  1 / 7  +M  ,  n^i-x)- ^  = (-n^x)- ^ =  1  1  and  F ( a + 6 x) = [1 - n - ^ - x ) - ^ - expi-i-x)- ^), n  1  n  x < 0.  1  n  For 7 = 0 , from Theorem 2.1.3 of Galambos (1987), we find a — logn, n  so that  F(x)  = e~  x  F(a + b x) = n ^ e * -  n  n  and  F ( a + b x) = [1 - n - V * ) " = e x p ( - e - ) . n  x  n  n  10  b = 1 n  The G E V distribution is useful for statistical inference when the distribution is unknown, and the tail index must be estimated. Since min(Xi,..., X ) - - max(-AV,..., - X ) , n  n  then, similarly for minima there are also three types of nondegenerate distributions as the univariate limiting extreme value distribution which can be summarized together and can be expressed in terms of (2-4-1)  L(x)  = =  l-H(-x) 1 - exp{-[max(l +  ~ ^ , 0)]- }. 1/7  <j  2.5 Multivariate Extreme Value Distributions We call a p-dimensional distribution function -P(x) nondegenerate if all of its univariate marginals are nondegenerate. One object of the multivariate extreme theory is to seek conditions on F(x), under which there are sequences {a }  and {b } of vectors such that each component of {b„}  n  n  is positive and < a + b z ) = H n ( a n + b „ z ) -»• H(z)  P(Z„  n  n  n  (2-5-1)  n  for a nondegenerate p-dimensional function H(z). Suppose X i , . . . , X  n  are independent random vectors, with common distribution  function F, then  H {z) n  = P(Z  L (w) = P(W n  n l  n l  < z . . . , Z < z) = u  n p  > ,...,W >w ) Wl  np  p  p  =  where G(w) = P(Xi > wu ..., Xp > wp). 11  F {z), n  G (w), n  Any problem on W is equivalent to one on Z by changing the basic vector X to (—X), therefore we concentrate on the vector Z of maxima. Next we list and discuss various properties which F(x.) or H(x.) may possess. P r o p e r t y 1. If F ( x ) n  is a sequence of p-dimensional distributions, let the «th uni-  variate marginal of F „ ( x ) be Fjf>(xi).  If F (x) converges to a nondegenerate continuous n  distribution function -F(x), then for each i with 1 < i < p, F ^ ( X J ) converges to the ith marginal F^(xi)  of F(x).  Later we shall see that all limiting distribution functions of multivariate extremes are continuous. Hence, Property 1 tells us that we can appeal to univariate case for determining the components of a and b whenever n  (2-5-1)  n  holds.  The following several properties are important in multivariate extreme value theory. P r o p e r t y 2. Let X i , . . . , X be iid p-dimensional vectors with common distribution n  function F(x.).  Then there are vectors a and b > 0 such that ( Z n  n  to nondegenerate distribution function H(x.),  n  —a  n  )/b  n  converges  if and only if, each marginal belongs to  the G E V family, and if C [y\'\...,y J ] n  l  F  n  -> CW(vi,• • •\y ),  n^oo.  P  Property 2 tells us that given a p-variate distribution function .F (x), we can check n  if its marginals belong to the G E V family, if so, then we use the methods of the univariate case to determine the components of the vectors a and b ; furthermore we can n  n  determine CV(y) by its definition and check if C ^ y / " ) converges. 1  Example 3. Let (X, Y ) have a bivariate exponential distribution  F ( x , y ) = l-e- -ex  y  + G(x,y),  where  G(x,y) = P ( X > x , Y > y ) . 12  (2-5-2)  If ( Z - a ) / b n  n  converges to a distribution H(zi,z2),  n  (logn,logn), and b  then we can choose a  n  = (1,1), and obtain  n  e~  + e~  Zl  Z2  F ( l o g n + 2i,logn + z ) = 1  h G(logn + zi,lo'gn +  2  z ). 2  n For  the Marshall-Olkin distribution  G(x,y)  = exp[-rc — y — Xmax(x,y)],  The  last term of (2-5-2) is  exp{  — (2 + A)logn — z — z — A max(zi, 2 )} x  1  =  2  2  A > 0.  ex'p{—z\ — z — \  2 + A  ri  =  ma,x(z ,z  2  1  1  0(n-"- ), A  and lim  F ( l o g n + zi,logn  + 2)  n  2  =  lim  n—>oo  e  1  - 2 i _(. -^2  \  e  TI  n—>oo V  = For  /  exp(—e * — e -  1  - 2 2  71  I  ).  \  another example, let  G{x,y)  = (e  6x  + e - l) 6y  _ 1 / e  ,  6> > 0.  Then G(log n +  log n + z2)  = (nV „  2 1  +nV  I( ^i n  2 2  l)~  -  1/e  )"!/«,  e  and lim  F (logn + ,logn + * ) n  2 l  2  =  lim(l-^—^  =  {exp[-(e + e )]}{exp[( ^ + e* )" /*]}. 21  -  K  ^—L  22  22  e  13  )» 1  2  In the second example, there is dependence in the limiting extremes. Property 3. A p-variate continuous distribution function H ( x ) is a limit distribution in (2-5-1) if, and only if, its univariate marginals belong to the G E V  family and if  its copula CH satisfies  C {y\'\ k  H  ..., yl )  = C (yi,  lk  H  • • •, yP), k = 1,2, • • •.  (2-5-3)  Property 3 tells us that if the limit H ( x ) exists, then we check condition (2-5-3). If it holds and if C'H is a copula, then we have got the actual limit distribution. Example 4.  Let H\(x),...,  H ( x ) belong to the GEV p  H(x)  =  family, then  H {x)..-H {x) x  p  is a possible limit of (2-5-1). This limit consists of independent univariate marginals. By assumption, the condition of Property 2 on the marginals is satisfied. Furthermore, by definition, CW(y) = y\ • • • yP for which the condition of (2-5-3) is evident. A n appeal to Property 3 yields the claim. Example 5. The distribution function H(x!,  ...,x ) p  = exp{exp[- m i n ( x i , . . . ,  x )}} p  is a limit in (2-5-1). The marginal distributions are Hi(x{)  = exp(—e~ ) = ^ ^ ( x , ) . Therefore it remains Xi  to check the validity of the condition (2-5-3): exp{— exp[— m i n ( x i , . . . , x )]} p  = min{exp[— exp(—x,-)]  and  C (yu H  • • •, y ) = m i n ( j / i , . . . , y ). P  p  14  1 < i < p},  Hence for k > 1,  C (yl , k  yl )  ,k  = [mm(yl ,  lk  H  vl )]  /k  /k  = m i n ( , . . . , y ).  k  p  y i  Example 6. The distribution  H ( x x ) = # ,o(si)#3,o(>2)[l + ^(1 - #3,of>i))(l - #3,o(z ))] u  2  3  2  does not occur as a limit in (2-5-1). (See Galambos(1987).) Even though the marginals  H\(xi)  =H  3 f i  ( x ) and  the copula  i / 2 ^ 2 ) = -#3,0(22)  1  2/2) = ym[l  C {yi, H  - 2/2)]  + ^ ( l - yi)(l  fails to satisfy (2-5-3). Let  X = (Xi,...  be a vector with distribution function F ( x ) . Let j(k)  iX ) p  =  (ji, • • • ,jk), 1 < k < p be a vector with components 1 < j i < J2 • • • < i t < P- The distribution function Fj(k)( h,  • • • ^j*) is  x  is obtained from F(x)  a  A;-dimensional marginal distribution, which  by letting X{ —> +oo for all % £ {ji, • • • ,jk}-  Let Gj{k){ JIT  ••i jk)  X  x  =  F  (-^ji  >  x  j n • • • •>Xj  k  >  Xj ). k  Assume that F ( x ) is such that each of its univariate marginals belongs to the G E V family. Then, there are a, and 6, > 0 such that n  F ( a , „ + binx)  lim where Hi(x)  n  n  t  = Hi(x),  1 < i < p,  (2-5-4)  belongs to the G E V family. We assume that a, and 6 n  determined and put a P r o p e r t y 4.  n  = ( a i „ , . . ., a ) and b pn  (Z — a )/b n  n  n  n  = (b ,..  tn  > 0 have been  ., b ).  ln  pn  converges to a nondegenerate distribution i / ( x ) if  and only if for each fixed vector j(k) and x for which Hi(xi),  1 < i < p, of (2-5-4) are  positive, the limit nHS> Gj(k){ jin n  a  ~i~ bj Xj , in  in  . . . , Oj n fc  "H bj n.Xj ) k  15  kn  — hj (Xj k  1  , . . . , £j ) t  (2  5  5)  are finite, and the function  H(x)  = exp{£(-l)  h (x ,... )}  £  f c  ih  h  (2-5-6)  lXjk  1 < J 1 < J 2 — <3k<P  fc=l  is a nondegenerate distribution function. If the actual limit distribution of ( Z — a ) / b n  n  n  is the one given in (2-5-6), then the following inequalities hold. Let s > 0 be an integer, then #(x;2s + l) < H ( x )  < H(x;2s),  (2-5-7)  where i7(x, 0) = 1 and #(x,r) = e x p { £ ( - l )  h (x ,...  £  f c  jk  When s = 0 in (2-5-7), the inequality H ( x ,  ,x )}.  h  (2-5-8)  Jk  1) < #(x) < 1 obtains. This means that  iJ(x) is never exceeded by the product of its univariate marginals, and H ( x ) exceeds 1. Example 7.  . Let F(x, y) = 1 — e~ — e~ + G(x, x  G{x,y)  = (e  y  6x  +e  y) where  -I)- ' .  ey  1 0  From the univariate case we know that  lim i*7(logra +  x) = exv(-e~ ),  x > 0,  x  lim F "(log n + t/) = exp(—e ), _2/  2  Since G\(x) = e  - x  and G ( y ) 2  = e  _ y  War) = lim nG(logra + x) = lim n—»oo n—*oo Similarly, h ( y ) 2  y > 0.  ne~^  n + x )  = e~ . x  = e , and _J/  h (x,y) 12  =  Jim nG(logn + x,logn + 16  y)  never  lim n(nee6x  =  +n e e  n—*oo n—+00  — l)~  1/,fl  JI<>  (e * +  e )- ' .  e  =  9 y  0y  1 e  Therefore, the limiting extreme value distribution is H(x,y)  = exp{-h (x)-h (y) 1  =  +  2  h (x,y)} 12  e x p { - e - - e~y + (e8x + x  e )~ }. By  1/e  The following definition is needed for the next property. Definition: The p-dimensional unit simplex S is the set of vectors q with nonnegative components qi such that YA=I Qi = 1Property 5 (The Pickands Representation for a min-stable exponential distribution). G(xi,...,  x ), with univariate exponential margins, is a survival function satisfyp  ing — \ogG(txi,...  ,tx ) p  = —t\ogG(x ,... 1  ,x ) p  for all t > 0 if and only if G has the representation -logG?(a;i,...,a;p) = /  [max. (q x )]dU(q ,... i  Js  p  where S  p  — {(qi,---,q ) P  i  1  ,q ), p  x, > 0,i = 1,... ,p,  (2-5-9)  1<»<P  ' qi > 0 , i = 1,... ,p, Yli qi — 1} is the p-dimensional unit  simplex and U is a finite measure on S . p  Property 5 is a result of Pickands and an alternative statement of Theorem 5.4.5 of Galambos (1987). It is not an easy task to give the representation (2-5-9) for a given H(x), but one usually does not aim at giving another form of H(x) when it is already known. The value of (2-5-9) lies in its possibility of generating functions which are limits in (2-5-1). For applications of Property 5, see Joe (1990a). A n example from there is:  exp{-A( , Zl  . . . , z ; \ , B c ( l , . . . , p ) , 6 ) } = e x p { - £ A ( £ x ) ' },6 6  p  s  B  t  l  s  > 1. (2-5-10)  17  For each integer m > 2, it can be directly shown that  x m ) = e x p { - ( ^ + --- + ^) / } 1  (2-5-11)  5  is a survival function over Xj > 0,j = 1,... ,m,for 1 < 8 < co. (2-5-11) satisfies the homogeneity of degree one condition of Property 5. Hence by Property 5 , if m < p, there is a measure V  m  on S  m  and a measure U  on S such that  m  p  where Um is a representation of Vm in p dimensions which puts all mass on qi,... ,q, Therefore, £ M £ ^ )  B^0  1  /  f  i  ieB  has the representation  where dU has the form  J2  B  \BdUs-  B y Property 5, (2-5-10) is a distribution.  2.6 Some relation between extreme value theory and Central Limit theory Since Central Limit theory is more familiar to statisticians, in this section we show some similarities and differences between Central Limit theory and Extreme Value theory. Consider the univariate case first. Suppose that X\,... ,Xn are iid with distribution  F. Let  Z  n  W = mm{Xi,..  = max{Xi,...,X },  n  n  18  .,X }. n  We  would like to know whether there are any sequences { a } n  S —a n  — A  Z- a  n  n  '  and {b } such that n  W — a 1  n  n  A  '. °  r  n  converge in distribution. From Central Limit theory, if F has a finite second moment, then there exist and {bn}  such that the limiting distribution (Sn — an)/bn  {a } n  has a normal distribution.  Definition. Let H ( x ) be a nondegenerate distribution function, then F(x) domain of attraction of H ( x ) if there are sequences { a }  is in the  and {b } > 0 such that  n  n  lim F ( a + b x) = H(x). n  n  n  From extreme value theory, we know that the limiting distribution H ( x ) of ( Z — n  o-n)/b  n  belongs to the GEV  attraction of For  family if the distribution function F(x)  H(x).  the multivariate case, assume that X i , . . . , X  X; = (Xn,..., Z Are  is in the domain of  n  n  are iid with distribution F . Let  S =Xi+••• + X ,  Xip),  n  = ( m a x ( X n , . . .,X ),... nl  n  ,max(A"i ,..  .,X )).  p  np  there any {a } and {b } > 0 such that n  n  converge in distribution? If the limit of (S — a ) / b and the second order moments of n  n  n  F exist, then from Central Limit theory the limit of (S — a ) / b has a multinormal n  distribution.  19  n  n  If the limit of ( Z — a ) / b exists, from extreme value theory, each univariate margin n  n  n  must be in the G E V family. A n analogy between multinormal and multivariate extreme value distribution is: (1)  Suppose ( Z i , . . . , Z ) has a multinormal distribution, then linear combinations p  of univariate normal, i.e. a\Z\ + • • • + a Z has a univariate normal distribution. p  p  (2) Suppose Xj ~ H(-,~fj, fj,j,o-j). After transformation to an exponential survival function Zj = -  ^og{H(Xf^ n aj)) h  h  have an exponential distribution. If (Z\,..., Z ) has a min-stable multivariate exponenp  tial distribution, then from Property 5, G(zi,...,z ) p  = exp{-A(z ,...,z )} 1  (2-6-1)  p  and V =. m i n ( ^ - , . . . , ^f-) has an exponential distribution - ( i>-> p) ^ where Wj > 0. vA  w  w  e  That is, weighted minima are univariate exponential. Note also the main difference: The multinormal distributions form a finite-dimensional parametric family, the multivariate extreme value distributions form an infinite-dimensional family (as given by Pickand's representation). Since parametric inference is easier than nonparametric inference, a "goal" is to find good parametric (finite dimensional) subfamilies of the infinite-dimensional family (see the next section). 2.7 Parametric families of multivariate extreme value distributions. A multivariate extreme value distribution with margins transformed to a survival function with exponential survival functions as univariate margin is a min-stable exponential distribution. Let G be a p-dimensional min-stable exponential distribution, from extreme value theory A = — log(G) satisfies A(tzi,...,tz ) p  = tA(zi,...,z ), p  20  Vt > 0.  Let G be the marginal survival function of G, then G ( z ) = e x p { — A ( z ) } where A S  s  a  s  s  S  is obtained from A by setting Zj — 0 for j $ s and s is a subset of {1,... ,p}. Here we will list some parametric families of multivariate min-stable exponential distributions in terms of A . These families can be found in Joe (1990a,b), Tawn (1989), Coles and Tawn (1991), Husler and Reiss (1989). The interpretations of the parameters are given after all the families are listed. The families are given here for p = 2,3,4, from which the general multidimensional form can be seen.  A(z z ;S) A( ,  s  6  1  r  2  x  2  r 2  )~  x l  l  3  (2-7-1)  r > 0,  \  + zi*) l \  A ( z z2, z3; 8, 8 ) = ((z* + zs2)s*'s u  1 s  2  z ; r) = z + z - (z7 + z  ZL  8>1,  = (z + z ) ' ,-  u 2  s  6 > 8  (2-7-2) (2-7-3)  > 1 ,  3  A(ZI,Z ,Z ;T,T ) 2  =  3  3  z + z + z - (z~ + z ) T  1  2  ( ( ^  + ^  T  T  )  T  2  3  1 / T  2  3  /  + ^3  T  A ( z z , z , z ; 8,8 ,8 ) u  T  3  4  3  T 3  )"  1 / T 3  - (z7  5  1  1  2  r>r >0,  ,  (2-7-4)  3  = ([(z + z ) *l 6  s  6  2  4  + z ^ ) - ^ - (z ^ + z ^ ) ' ^ +  T3  + zi*} *' > + zi*) ' *, s  s  1 8  l < 8  4  < 8  3  < 8 ,  (2-7-5)  A ( z , z , z , z ; 8 , 8 , 8 ) = {{z[ + zs2)Si'& x  2  3  2  4  4  + (z > + z ^ * ) ' * , s  3  6  1  6  1<  6 ,8 2  4  < 8,  (2-7-6)  A(zi,z ,z , 2  3  z; 4  T,r ,T ) 3  4  21  =  z + z + z + z - (z7 + z r  1  2  (z7  3  +z7 )-  Ti  T i  3/T  {{ r +*3  r  T  T3  Z  (((*r + ^ ) T  {z *  -  1 / T i  (c*r + ^ T 3  T  4  T 4  4  +z " ) T3  T3  )"  1 / T 4  T 3  ~  +^"  1 / T 3  T 4  )  -(z  T  + ((*r +zrr  1/T4  3  + ^" )"  T4/T3  T4  3  1/T4  ,  4  r  i/T3  +z7 )~  T 3  T 3  2  +  1/T4  +*r )~ 4  1/T4  r > r > r > 0. 2  -  1 / T 3  +  _ 1 / T 4  + (c*r+*-2 r' +*r )~  1/T3  +*r )~  - (z~ + z ^ ) -  1 / T  4  T3  T 3 / T  )~  +2 "  T  2  +^- )"  /T3  2  (2-7-7)  4  A(ZI,Z ,Z ,Z ;T,T ,T ) 2  3  4  2  4  = z + z + z + z - (z7 + 2 - T T  1  2  3  (zp + z r ) 4  (c*r + * n 4  -T2  3  O r  +* "T 2  - (sp + z - T T  T4  1/T4  2  +*3- )" T 4  + ^r ) 2  4 / T  T4/T2  1 / T 4  -  +* r r  + ((*r +* - ) T  T 4 / T  2  )"  + (*r + (*  1/T4  4  T 2  T 2  3  _T2  3  + (*" + - )  - ( ^ +^  1/T4  3  - (^ " +*rT  _ 1 / T 4  T 4 / T  (*r + (*  1 / T  2  4  T 4 / T 2  2 4  )  _ 1 / T 4  1 / r i  +*r )~ 4  + *r ) 2  2  4  1 / T 4  +  )~  1 < r , r < r.  ,  )-  +  1 / T 4  T4/T2  T 4  1/T4  -  (2-7-8)  Finally, the Hiisler and Reiss(1989) model with p = 2 is A(zx,  A) = *(A +  l  0  g  Z  l  -2A l  0  g  Z  2  ),  + *(A + ° 1  2  g Z 2  -2A  l 0 g  ^), , 1  A > 0. (2 - 7 -  For p > 2, let Xij € (0,oo)  VI < i , j < p with  i ^  Put  Moreover , for 2 < k < p and m = ( m . . . , m^) with 1 < m i < m < . . . < 1 ?  2  <  define  rfc,m —  ^ ( A ^ . ^ + x . — \ .^ .) i,j<k-\ 2  2  2  m  mk  m  22  m  (2-7-H  Furthermore , let  denote the survivor function of a normal random vector with  S'(-|r)  mean vector 0 and covariance matrix T. The extension of (2-7-4) is  A(z ,...,z ;\ j, 1  p  i  \ < i < j <p)=  p  k  J2(~ ) '52 k,m(-1ogz ,...,-\ogz ) 1 k  h  mi  k=\  mk  (2-7-11)  with ,  /•oo ^m(y)  =  /  S  -,  (y- - Z + 2 A f  2  jf- |r ,m 1  n > ) m  f c  6  ~^>  Vk  J  for 2 < k < p, where J2 nieans summation over all p-vectors m = ( m i , . . . , m^) with k  1 < mi < m < ... < m 2  < p, and /ii, (y) = e~ for m = 1,... ,p. The cases p = 3,4 y  m  k  are given explicitly in Chapter 3. For comparison, note that the bivariate models (2-7-1), (2-7-2) and (2-7-9) each have a single dependence parameter. For (2-7-1) and (2-7-2) dependence increases as the parameter increases; for (2-7-9) dependence decreases as the parameter increases. The trivariate models (2-7-3) and (2-7-4) have respectively models (2-7-1) and (2-7-2) for each of the three bivariate margins. For the trivariate model (2-7-3) (respectively (2-7-4)), 8 (respectively T) is the dependence parameter for the (1,2) bivariate margin; £3 (respectively, T3) is the dependence parameter for both the (1,3) and (2,3) bivariate margins. We have something similar for the 4-variate models (2-7-5) to (2-7-8). The 4-variate models (2-7-5) to (2-7-8) also have respectively models (2-7-1) and (2-7-2) for each of the six bivariate margins. For (2-7-5), 8 is the dependence parameters for the (1,2) bivariate margin; <$ is the dependence parameter for the (1,3) and (2,3) bivariate 3  margins; 84 is the dependence parameter for the (1,4), (2,4) and (3,4) margins. For (27-6), 8 is the dependence parameter for the (1,2) bivariate margin; 82 is the dependence parameter for the (3,4) bivariate margin; 84 is the dependence parameter for the (1,3), (1,4), (2,3) and (2,4) margins. The model (2-7-11) has a dependence parameter for each of the p(p— l ) / 2 bivariate margins. For (2-7-11), the dependence parameter of the 23  bivariate margin, with i < j, is A,j, and the  margin has form (2-7-9).  Models (2-7-3) to (2-7-8) and their extensions have p — 1 parameters in total; as above, some of the bivariate margins have the same dependence parameter.  Hence  these models do not have as much flexibilty in the dependence pattern but they have a simpler form than (2-7-11). 2.8 Estimation for a parametric family of multivariate extreme value distribution. Joe (1990b) suggests the following procedures for fitting of a parametric multivariate distribution to iid X i , . . . ,  X . n  (1) Fit the p univariate margins separately by maximum likelihood using the G E V family. (2) Transform each margin so that each transformed variable has exponential survival function Gi(zi)  = e~ . Zi  (3) Suppose that the transformed data are iid p-vectors, then compare fit of different families of min-stable multivariate exponential distributions. Check for parameter estimation consistency with bivariate and higher order margins. (4) Go back to the original data and estimate both multivariate parameters and parameters of the univariate margins simultaneously, using maximum likelihood if a copula has been decided on as being acceptable. In  (3), to check the fit of a model, we compare parametric estimates of A with  a nonparametric estimate of A. if (Zi,...,  Z) p  Besides the consistency checking from we know that  is a random vector with the min-stable exponential survival function  G — e~ , then A  E(wm[^,...,^]) Wi  for  =  w  [A{w ,;..,w )]-  1  1  p  p  u>j > 0, j = 1,... , p. So if Y i , . . . , Y  p  24  are a transformed random sample, expo-  nential probability plots of Vi = minjYij/wj  can be used to check the assumption of  min-stability for multivariate exponential distribution. If the plots are adequate, a nonparametric estimate of A(wi,...,  w ) is given by n/ p  parametric estimates of A from different models.  25  J27=i ^  which can be compared with  Chapter 3 Maximum Likelihood Estimation With The Husler-Reiss Model Hiisler and Reiss (1989) obtained a new class of multivariate extreme value distributions by taking a non-standard extreme-value limit. Let X and M  n  be as defined in Chapter 2. For a mutivariate normal random vec-  tor, having all correlation coefficients smaller than one, Sibuya (1960) proved that the marginal maxima M i , . . . , M „ n >  iP  are asymptotically independent.  The rate of con-  vergence to asymptotic independence is slower as the correlation coefficients increase. Since in practice the sample size n never goes to infinity but is a fixed positve integer, is there another asymptotic formulation which may provide a better approximation? This motivated the work by Hiisler and Reiss. For the bivariate case, Hiisler and Reiss let the correlation coefficient p — p increase n  towards 1 as the sample size increases. It is shown that the marginal maxima are neither asymptotically independent nor completely dependent if (1 — p )\ogn n  converges to a  positive constant as n —> oo. This was extended to dimensions p > 2. Hiisler and Reiss define their limiting distribution by the c.d.f. or survival function. To get a maximum likelihood estimation, the density is needed. To obtain the density, the following result for the conditional distribution of a multidimensional normal distribution is needed. Result 1: Let U = (Ui,..., U ) and V = (Vi,..., H ) be random vectors such that a  MVN  12  0, V  ^ £21  ^22  \  JI  where E n , £12, £21, £22 are a x a, a x b,b x a,b x b matrices respectively, and suppose £ 1 1 and £22 are nonsingular. Then U | V = v ~ MVN  ( EjaE-iy; E 26  n  - EiaE^i )•  Let F x and F x denote the distribution function and survival function for a random vector X , and let fx and /x|Y denote the density function of X and the conditional density function of X when Y is fixed. Then the joint distribution of X and Y is /x,y( »y) = /x|Y(x|y)/v(y)x  Let <^(x;E), $(x; E ) denote the multivariate normal density and survival function with zero mean vector and covariance matrix E . Then /u,v(u,v) = ^(u -  S H S ^ V J S H - S12E- E21)^(v;S22) 1  and  Fu,v(u,v) = $((u,v);E)  J  roo f  "1 '  /-co roo •••/  /  Jvi, Jui  too  •••/  0(s-Ei S - t;Eii-S Sj E )^(t;E )dsdt 1  2  Ju  2  1  2  1 2  2  2 1  2 2  a  •••/ ^ ( u - E i a E ^ t j E u - E w E ^ E a O ^ E M j d t .  Therefore, 3 <&((u,v);£) _ dv\ • • • dvb b  (  _  1  )  6  $  (  u  _  S  i  2  S  - i  v  ;  E  i  i  _  E  i  2  E  - i  S  2  i  )  ^  (  v  ;  ^  Now we are ready to derive the density for the Hiisler and Reiss survival function when the univariate margins exponential have been transformed to distributions with mean 1. We provide details for dimensions p = 2,3,4 from which the general case will be apparent. The notation gets increasingly difficult as p increases. 3.1 B i v a r i a t e Case Frome Section 2.7, G(wi,w ,X) 2  = e x p { - A ( i 0 i , u ; , A)}, 2  27  2  where  A) = $(A + ^ \og(wi/w ))w  A (w ,w , 2  1  2  2  Since J^- = — f ^ - G and J^- = — dG 2  =  _d_  dw\dw  2  a  =  1  2  x2  -\og[cf)(\+  dA  dw  2  From 4>(x) = (2Tr)~ ^ e~ ^  A > 0.  2  we obtained the density function of t7(u>i, w ,  dwi dw  2  + $(A + ^ log(w /u;i))u)i,  2  1  2  &4 g A 2  =  dw  2 c  dw\ dw  2  _  d A 2  2  dw dw  2  x  2  we have ^'(x) = a;^(x). Also  2  ^\og(w /w ))/w ] 2  1  1  =  A / 2 + (logto - logwi)/2 + ^ - ( ( l o g w - l o g ^ i ) -t-logiu! + log V^TT  =  A / 2 + (log w  =  - l o g t y ( A + ^-log(tOi/u; ))/u; ].  2  2  2  2  2  - log w )/2 + g^-((log tux - log w )  2  x  2  2  2  + log w  2  + log V2TT  2  Hence <KA + ^ log(u; /u)i))/u;i = <^(A + ^  log^/u^))/^.  2  Therefore  dA  2  d 1 1 d w ^ ($(A + 2A — l o g ^ '/ u ^") ) ^ + $(A +2A — log^/w^w-^) v  =  ~ ^ (  A  toV  + 2Xlog(w2M))(u> M) + ^ ^ (  +$(A + =  v  ^(X + 2A  2  ^\og( /w )) Wl  2  ^-\og( /w )). Wl  2  Similarly,  3\A.  1  Finally,  28  A  + ^-log(^iM))(^i/^2)  d o^  m X + 2  ^ ( A From Ait2,  1 2X  l o g { w i / w 2 ) )  + ^logKM)).  A ,i and A\ we can get the density function g(w ,w , 2  2  1  A).  2  3.2 T r i v a r i a t e Case Let f  E =2  2 \ \2 _i_ \2 ^13 ' ^23 "12 2A^3 /  9\ 13 \2 I \2 \2 13 "t" ^23 12 2  y  A  — A  be the matrix in (2-7-10) and let \2 _ ^13 + ^23 ~~ \2 A  P12 = P21 -  P13  = psi =  2A13A23  1  \2 , \2 12 + 23 23 1 2Aj A23 \2 . \2 12 + 13 23  A  A  A  A  A  2  P23 = P32 =  —  A  2A12A13  i.e., Pij =  2A,/;Ajfc  for i,j, k distinct. Then h 3 < 1 2 3 in (2-7-11) can be written as  C(wi,w ,w , 2  3  A12, A 1 3 , A 2 3 )  = / JO  $ ( A + — — ( l o g g - l o g t u i ) , A 23+777—(log g-log w2); LA ^ 13 13  A  23  where with some abuse of notation 1 1 P  1  $(•;/>)=¥(•;  \  P 1 From Section 2.7, we have  G(w ,w ,w ) 1  2  3  - exp{-A (w ,w ,w , 3  1  2  29  3  X , X , A )} 12  13  23  p )dq 12  where with simplification of (2-7-11),  A (w ,w ,w ,Xi ,Xi ,\ ) 3  1  2  3  2  3  ~w\ - w - w + A (wx,w ,  23  2  3  1  Ai ) +  2  2  A13) + A (w ,w ,  A (w ,w , 2  2  3  2  2  X )+  3  C(wi,w ,w ).  23  2  3  Constraints are 0 < A < 0 0 and 0 < pij < 1. tJ  By symmetry, equivalent forms for C can be obtained by permutations of subscript indices: C(iui, w , w , A , A i , A ) 2  3  rw\  = /  12  3  23  1 ^(Ais + 77T-(log q - log w ), A 3  23  1  $(^12 + 7 T T ( S 9 - logioi), A -  0  2  1  ru>2  J'  1 + -r— (log q - log w ); p )dq  12  + — — (log 9 - log w );p )dq.  lo  6Ai2  The density function for G(wi,w ,w ) 2  2 3  3  is:  3  13  4A23  d  3  g(w ,w ,w ) 1  2  3  = ( - l ) ^ — 7 ; — « — e x p { - A ( w , w , w , A i , A i , A )} 3  3  OWiOW OW _ dA d A dw\dw dw 2  _ From  SMdMdM dwidw dw 2  2  1  2  3  2  3  23  3  _ 8A d A dw dw\dw  2  3  2  3  2  3  3  3  2  3  _ dA d A dw dw\dw  dA dw dw dw  2  3  3  3  3  +  3  2  1  <j>(xi,x ; p) = (f>( )(f)(^=^) 2  Xl  00  f°° ±t i ~ P ~L z  Z  1  and ^ , x ; p ) = 1  2  $ (  7  f  = ) ^ i ) .  Let Cl,23 :=  = $(Ai3 + atvi  \og(Wi/w ), A12 + J-J—log(lWi/u7 );P23)zA 3  zAi3  2  1 2  30  2  3  B y symmetry, we have C ,i3 •- dC and C ,i2 := 2  Similarly,  3  dC 2  12,3  dw\dw<i  $[(( i3 + 2^1og(u7i/u;3)) A  - ^ (Ai2 + ^ M t o i / w j ) ) ) / ^ 12  - P23}  2 3  x  <K i2 + ^—log(u;i/u;2)); A  2Ai2io  B y symmetry, we also have C  2  := g^f^ and C , i := afffe"- Finally, 23  1 3 i 2  dC 3  C i o :— 3  dwidw dwdw 2  3  3  1  <K 13 + ^T—log(w M),A 3 + A  3  £A\  ^T—^Og(w /w2);pl )jT-^ 6A23  2  3  3  <K i3 + r r  2  ^^13  Combining C  i  ±A\ A2 W\ 3  3  t«2  log(^3/wi))<^[((^23 + 7 7 — l o g ( w / w ) ) -  A  3  ^^23  13 log(wi/w3)))/y/l + ^77  Pi2(Ai3  2  - p\ ] X 2  4A  A  1 3  ^2 \ / l ~ ^ - P i 2  1 23  , C ,i3, C7,i2, Ci , , Ci , , C ,i, Ci23,and f^ -, f ^ f ^ ^ g ^ , 2  li23  2  3  9~~9^7> and QJ^Q^£Q  w  2  3  3  2  23  a£rfc  have an expression for the density function <7(it>i, iw> u> ).  e  2  W  3  3.3 Four-Dimensional Case Let _ \2 + A2 24 12  \2 1 \2 _ \2 \2 1 \2 _ \ 2 14 + 24 12 14 + 34 13 \2 1 \2 _ \2 2A 24 + 34 23  -r1 A\2  \2 A  2A? E = 2  A  y  2 4  \2 A 14  A  4  A  34  \2 —A 13  A  A  A  A  24  24  A  1A \2 \2 —A 34  23  A  A  2A|  A  4  be the matrix in (2-7-10) and let Pijk  A  \2 1 \2 _ \2 t j + ifc jk ij ik A  A  ZA  A  where A,j = A - and, «, j , A; are distinct and between 1 and 4. Then / i i Jt  4 )  can be written as:  31  2 3 4  in (2-7-11)  D(w ,W ,W ,W ) 1  2  3  =  4  rm l / $( 14 + 7T\ ^g(q/w ),  l A + p-r—  A  24  4  JO  &A\  A34 + 7 7 7 —  ZA34  log(q/w ), 4  ^^24  4  \og(qjw ),  R(p i2,  4  P423))dq,  0413,  4  where  f 1 a 0A R(a, b, c) —  X  a l e 6 c 1  By symmetry, equivalent forms for D are  D(lWl,l0 ,f>3,tU4) 2  y^i  =  /  JO  1 $(Ai2 + 7 7 7 — log(?/w ), A 2  A^Vl  A14  +  rw2  2A 14  Jo  1 + 7 7 7 — log(o/w ), 3  ZA13  log(a/to ), R(p 3, 4  0124,  12  1 *(Ai2 + 7 7 7 — log(o/wi),A  /  a 3  23  P\3A))dq  1 + 777—log(?/w ), 3  £>A\2  ^ 23 A  •^—\og(q/w ),R(p i3,p2i4,p234))dq  A4 + 2  2  4  ZA2  4  /  ru>  Jo  3  \  \ •  $(Ai3 + 7 7 7 — log(g/wi), A + ^ 7 — log(?/u;2), AA\ ZA23 A34 + 77T— I°g(?/W4), #(0312, 0314, P324))dq. 2 3  3  ZA34  Constraints are: 0 < Xij < 0 0 , i,j,  i ^ j, 0 < /o.-jjk < 1,  are distinct and R(P412, 0413, 0423), R(P\23i 0124, 0134), #(0213, 0214, 023 ), #(0312, 0314, 0324) 4  are positive definite.  32  Taking derivatives, we have  dD  -Dl,234  —..  o^r*  1 ,  2^  { X i 2 +  i  ,  .  '  i  ^g(w /w ),  TT—log(«;i/u;3),Ai4 + 7 7 7 —  A13 +  .  l o g ( W i / W 2 )  1  R(p ,  4  123  0124,0134))  ZA13  By interchange of 1 and 2, 1 and 3, 1 and 4 in  1)1,234  we can obtain  As,i24 : = | £ , -04,123 := f j - B y making use of Result 1 with E 3u>3 :  /  1)2,134  1  := g ~ 0134  n  ^ 0134  1 j  \  0123  J12  , E22 = 1 we have  \ 0124 j  d D dwidw 2  D  12,34 • —  =  2  dD 234 dw lt  2  — 1 $( 13 + 7 7 T — A  ZA13  log(wi/w ) 3  1  ~ 0123(Al2  1  + 7 - — \ 0 g ( w 1 / w 2 ) ) , \ l 4 + 777—log(tt>l/w 4 ) ^Ai2  ^  1  i 4  1 1 ^°&( l /l)) 7y[  012 (Al2 + 7 7 7 — log(wi/w2)),#34,12)<KAl2 + 77\  W  4  ZA12  A  .  W  ^Ai2  ZA12W2  where 1  R.34,12  0134  ^ 0134 (  I  \  0123  0123 0124  1 j  1  ^ 0134  \  0123  —  —  0134  01230124  —  1  01230124  \  2 - 0124  )  By permutation of indices we can easily get  D 13,24 •— Q g. W  9D  9D 2  23  9wi dun '-^23,14  dw dw i 2  3  r) 24,13 •—  iy  2  dw^dvn  '  Again using Result 1 with E n — 1, / S12 —  \  0123  /  1  0134  0134  1  , E22 —  0124 )  I  33  \  /  dD 2  9 ^ 9 ^  5 -^34,12  -  5  then y-i ^22 — -#123,4 = 1 — E i S 2  1  1 1 -  Pl34  ^ —Pl34  E i = 1 — (p  2 2  -/>134  2  123  +p  1 2 4  1  j  — 2/)i 3/)i /) 34)/(l — p 2  24  1  1 3 4  )  and  d D 3  D 123,4  dw\dw dw 2  3  ${(Al4 + -^—logiwxIWi) - [(/9123 - />124Pl34)(Al2 + T"J— l o g ^ i / ^ ) ) + Z/\i4 ^^12 f>124 - /9l23/»134)(Al3 + ^ 7 log(li>i/u> ))]/(l ~  p\ ))I\J ™A) R  3  34  1 </>(Ai2 + -J—log(t«i/ti;2),Ai3 + r-J— l o g ^ / ^ a ) ; ^123)7. . ZA12 ZA13 4Ai2Ai3U;2^3 By symmetry and permution of indices we can easily get dwidufySwi D'  £*i24,3  Finally,  134,2 • dw^dw^dun1 -^234,1  d D A  D 1234  dwidw dwzdw 2  4  •^—log(w /w ),X  <^(Ai4 +  1  4  + -^—\og(w /w ),X  24  2  ZA\  4  + -^—\og(w /w ),) 3  34  ZA24  4  4  ^A 4 3  R{,P412, / > 4 1 , / > 4 2 ) ) g ^ ^ ^ ^ ^ ^ 3  3  From Sections 2.7 we have  G(w!,w ,w , 2  w ) = exp[—A (iwi, w , w , w , X , X , A , A , A , A , A )],  3  4  2  4  3  4  12  13  i 4  2 3  23  2 4  34  where  w , A , A , A , A 3, A 3, A , A ) =  A (w ,w ,w , 4  1  2  3  12  4  13  14  2  2  24  34  ^i(iwi) + h (w ) + hi(w ) +fei(iw )- h {w ,w ) - h^W!, w ) - h (w ,w ) 1  2  3  4  2  1  2  3  2  x  h (w ,w ) - h (w ,w ) - h (w ,w ) + h (u)!,w ,w ) + h (w!,w ,w ) + 2  2  3  2  2  4  2  3  4  3  2  h (u)!,w ,w ) + h (w ,w ,w ) + h (w!,w ,w ,w ) 3  3  4  3  2  3  4  4  2  34  3  4  3  3  2  4  4  =  IWi + w + w + w - [wi +w - A (w ,w , 2  3  4  2  2  1  A )] 12  2  [W!,W - A ( l U i , to , Ai )] - [wi + w - A (wi, 3  [w + w 2  3  2  3  2  2  3  23  [w + w - A (w ,w , 4  2  3  2  4  3  4  2  2  34  2  3  24  3  4  2  1  2  13  X , A ,A 12  23  14  A , A , A ) - D(u>i, 1 ^ 2 , ^ 3 , ^ 4 )  4  23  = -2(wi + u> + w + w ) + A (w ,w ,\ )] A (w!,w ,  4  2  C(w ,w ,w ,  i3  2  2  A )] + C(ioi,tw ,iW3, Ai , A , A ) + C(wi,w2,w4,  4  C(wi,ti> ,t«4,A ,A14, A 3 4 ) +  2  W , A14)] -  2  A (w , w , A )] - [w + w - A (w , w , A )] -  3  3  4  24  34  + A (wi,w ,\ )]  12  2  3  +  13  Ai )] + A (w , iw , A )] + A (u> , u>, A )] +  4  4  2  2  3  23  2  2  4  24  A (u; ,u;4, A )] + C(w , w , w , A , A 1 3 , A ) + C(u>i, w , w , A i , A , A ) + 2  3  34  1  2  23  C(w ,w ,w ,  C(iwi,io ,u;4, A 1 3 , A i 4 , A ) + 3  12  3  2  3 4  3  2  4  2  i 4  24  A , A 4 , A 3 4 ) - D(u>i, w2, w3, w4)  4  2  23  and g(w , ,w ) uW2 W3  =  4  ( - l )  aA aA <9A d A 4  4  4  ^  4 F  A  a A  i  8  ^ ^ ^ 3 , ^ )  4  aA4a^4  2  4  >  4  dw dw dw dwidw dw dw dA dA dA dA 8A dA dA 8A dA dw\dw dw dw dw\dw dw dw dw dw dw-^ dw dA 8A dA dA 8A dA dA dA dw dw dw\ dw dw dw dw\ dw dw\dw dw dw dA dA dA dA dA dA dw\dw dw dw dw\dw dw dw dw\dw dw dw dA dA dA dA dA dA dwidw dw dw dw dw dw dw dw dw dw dw dA dw\dw dw dw 2  3  4  2  2  3  4  2  4  4  3  2  4  2  4  4  4  4  2  4  2  4  3  2  2  4  2  4  4  3  2  3  4  2  4  4  2  2  4  3  2  3  4  4  4  3  3  |  4  4  3  4  3  4  2  4  3  4  4  3  4  2  |  2  4  2  4  4  2  4  3  4  4  3  2  4  4  2  4  1  3  4  4  4  2  2  3  4  4  x  A  4  2  3  4  Upon differentiation, dA  4  -2 + A  1>2  + A  1<3  + A  1A  + Ci,23 + Ci,24 + Ci,34 - D i ,234-  By symmetry,!^ ,1^-, f^ are like f^- with indices 1 and 2, 1 and 3, 1 and 4 inter1  1  changed. Also, Pi A. 1  dw dw 1  2  A\2 + Cl2,3 + Cl2,4  35  —  ^12,34  24  and o  Z  — W23 ~ -L/123,4  5  OW1OW2OW3 By the same reason as the above we can easily get dA  dA  2  dA  2  4  dw\dw  dwxdw  3  dA dwidw dw '  3  dw dw  4  3  4  z  4  4  4  dA dw dw dw  3  4  2  4  3  dA dwidw dw '  3  dA  2  4  dw2dw dw2dw  4  2  dA  2  4  4  4  2  3  4  .Finally, dA  ,  4  4  Q Q Q Q C'Zi;iC'l«2C'W3C'^'4  —  ~~ ^1234-  ;  3.4 Computer implementation of the Hiisler-Reiss model The parameter A,j's are estimated with a quasi-Newton routine. The multivariate normal survival function is computed with the routine of Schervish (1984). Other routines that are needed are a numerical integration routines and routines for the functions:  C i 3 , C i 2 , , Ci23,-Di,234, ^12,34,-Di23,4,-Di234i 2  3  Then C 1 < 2 4 , C 1 3 i 4 etc., can be com-  puted from these routines by changing the parameters of the routines. 3.5 Data Analysis Joe (1990b) constructed some multivariate extreme value environmental data sets from Bay Area and from the Great Vancouver Regional District. The first data set consists of weekly maxima of ozone concentrations ( in parts per hundred million) for several monitoring stations in the Bay Area; the weeks were for the months of April to October in 1983-1987. For comparison with Joe (1990b), the same subset of 5 stations are used. The second data set consists of weekly maxima of SO2, NO2 and ozone concentrations (in parts per hundred milion ) for several monitoring stations in the Greater Vancouver Regional District: the weeks were for the months of April to September to October for July 1984 to October 1987. Joe (1990b) eliminated all weeks with too many 36  missing values, and for comparison the same reduced data are used. In Chapter 4, weeks with missing values are used as censored data. Details of exploratory data analysis of both data sets are given in Joe (1990b). Joe (1990b) compares the fit of the models (2-7-1) to (2-7-8) and their multivariate extensions, as well as other models in Joe (1990a). The models in Joe (1990a) which have some appealing theoretical properties were found to fit worse even though the models had more parameters. A partial explanation is that the parameters in these models are not easily interpretable as in (2-7-1) to (2-7-11). For comparison with a baseline, the likelihood of the multivariate normal copula with exponential margins was computed. The comparison of the Husler-Reiss model with the best fitting of (2-7-1) to (2-7-8) and with the multivariate normal copula is given in Tables 1 to 7. Only the parameter estimates for the Hiisler-Reiss model are given. Note that the likelihoods given are for the dependence parameters, after the univariate margins have been transformed to exponential survival functions (see Section 2.8). Also all likelihoods based on margins ( or subsets of data from subsets of stations) are given. This is to check on consistency of the parameter estimates for multivariate margins of different orders. This check provides an indication of the goodness of fit of the multivariate models. A final note for the tables is that for 3 or more stations, labelled as l , . . . , p , the parameter estimates are given in lexicographical order, that is,  A12, A 1 3 , . . . , A i , A 2 3 , . . . , A p  2 p  , . . . , A _i . p  i P  For example,  consider the pairs of stations P T , SJ in Table 1, where the other stations are abbreviated C C , V A , S T . In the trivariate subset, ( C C , P T , SJ), the A parameter for (PT,SJ) is the third value (0.588); in the trivariate subset, (PT, SJ, V A ) , it is the first value (0.595). in the 4-variate subset, ( C C , P T , SJ, V A ) , it is the fourth value (0.567), etc. Recall from Chapter 2, that a smaller value of A means more dependence. Some overall conclusions from the tables are:  37  (1) The parameter estimates for a particular A are similar over different likelihoods. (2) Of the 37 bivariate likelihoods, the likelihood of the Hiisler-Reiss model is higher than those of both models (2-7-1) and (2-7-2) in 13 cases, and it is higher than that of the multivariate copula in 18 cases. (3) Of the 25 trivariate likelihoods, the likelihood of the Hiisler-Reiss model is higher than those of both models (2-7-3) and (2-7-4) in 11 cases. With the adjustment for the number of parameters ((2-7-3) and (2-7-4) have 2 parameters, the Hiisler-Reiss trivariate model has 3 parameters), through a penalty of half the number of parameters for the log likelihood( the Akaike information criterion), the Hiisler-Reiss model is only better for the first data set in Table 1. (4) Of the 9 four-variable likelihoods, the likelihood of the Hiisler-Reiss model is higher than those of the models (2-7-5) to (2-7-8) in 6 cases. This is true even with adjustment for the number of parameters ( 6 for the Hiisler-Reiss model and 3 for models (2-7-5) to (2-7-8) ). Overall, it appears that the Hiisler-Reiss model provides a better fit for Tables 1 and 3, and maybe Table 5. Therefore it is a useful model in addition to those in Joe (1990b). A n indication of whether the Hiisler-Reiss model is better can be seen from the estimates for the bivariate likelihoods. If the pattern of the A's does not fit the bivariate patterns that are possible with models (2-7-3) to (2-7-8), then the Hiisler-Reiss model should be better. ( The restrictions of the bivariate patterns for models (2-7-3) to (2-7-8) can be seen from taking the bivariate margins in these models.)  38  Table 1. Bay Area - Ozone(n=120) pay A Notation: Concord(CC),Pittsburg(PT),Vallejo(VA), San Jose(SJ) Santa RosT). Tikelihood for dependence parameters negative log likeli.hood estimate r^s subset T F T T normal H-R model CC, PT — O S 3 — 154.75 154.28 154.08 153.50 175.96 179.40 176.91 0.451 176.17 cd, SJ 190.30 0.532 190.70 190.38 192.41 CC, VA 210.47 211.60 211.03 214.24 0.708 CC, ST 199.06 199.85 199.34 200.48 0.595 PT, SJ 197.15 196.07 196.05 198.52 0.581 PT, VA 213.92 214.80 214.37 214.54 0.749 PT, ST 197.15 197.68 197.40 202.49 0.579 SJ, VA 213.52 214.56 213.96 218.74 0.740 SJ, ST 199.54 199.86 199.62 203.12 0.599 VA. ST 212.23 221.32 216.41 212.94 0.354,0.453 CCJ^SJ 0.5^88 222.92 CC,PT,VA 0.353,0.532 225.26 221.36 223.87 0.&0 243.86 CC,PT,ST 0.352,0.704 244.94 245.20 245.70 0.750 242.80 CC,SJ,VA 0.451,0.536 245.38 243.37 248.67 0.576 0.450,0.716 268.11 265.32 272.20 265.69 CC,SJ,ST 0.741 268.28 CC,VA,ST 0.533.0.708 271.62 271.40 273.54 0.606 0.595,0.583 267.09 265.27 269.90 264.27 PT,SJ,VA 0.579 287.04 PT,SJ,ST 0.594,0.755 290.80 288.08 291.56 0.?38 275.45 PT,VA,ST 0.582.0.752 279.04 277.97 278.94 0.d06 0.580.0.747 282.20 278.01 284.13 275.85 SJ,VA,ST 0.607 CC,PT,SJ 0.357,0.441 288.42 280.47 280.01 275.72 VA 0.539,0.567 0.581,0.574 CC,PT,SJ 0.355,0.441 311.54 305.04 303.61 299.69 0.70. 0.567 ST 0.738,0.741 CC,PT,VA 0.353,0.534 307.17 303.67 304.20 300.66 0.706,0.581 ST 0.754,0.607 CC,SJ,VA 0.450,0.538 329.94 324.70 329.57 321.01 0.714,0.577 ST 0.750,0.609 0.594,0.585 352.88 347.34 349.94 PT,SJ,VA 342.36 0.755,0.581 ST 0.746,0.610  39  Table 2. GVRD-Sulphur Dioxide(n=70) Notation: Confederation Park(CP), Second Narrows(SN), Anmore(AN),Rocky Point Park(RP). Likelihood for dependence parameters. estimate  negative log likelihood  subset  H-R  6's  T'S  normal  H-R model  CP, SN  0.845  129.55  129.54  129.37  129.58  CP, A N  1.086  136.11  135.96  134.82  135.94  CP, R P  0.905  131.45  131.67  129.61  132.00  SN, A N  1.068  135.36  135.54  133.97  135.62  SN, R P  1.229  137.04  137.40  136.48  137.70  AN, RP  0.756  123.39  124.54  123.63  126.08  CP,SN,AN  0.840,1.068  193.83  193.82  192.16  193.44  190.28  192.97  189.26  190.89  183.89  188.01  182.94  187.05  187.61  190.59  187.30  190.93  241.86  247.15  240.07  244.05  1.062 CP,SN,RP  0.846,0.891 1.157  CP,AN,RP  1.045,0.909 0.750  SN,AN,RP  1.064,1.172 0.750  CP,SN,AN  0.840,1.060  RP  0.900,1.032 1.155,0.748  40  Table 3. GVRD-Nitrogen Dioxide(n=84) Notation: Kensington Park(KP), Confederation Park(CP), Anmore(AN), Rocky Point Park(RP). Likelihood for dependence parameters. negative log likelihood  estimate subset  H-R  6's  r's  normal  H-R model  KP, CP  0.378  112.97  113.10  120.15  114.06  KP, A N  0.415  119.26  119.31  124.01  119.59  KP, RP  0.362  110.64 .110.64  119.83  111.53  CP, A N  0.466  127.84  127.30  128.81  126.55  CP, R P  0.346  108.53  108.27  107.87  107.71  AN, RP  0.356  110.48  110.22  112.64  109.91  KP,CP,AN  0.379,0.417  147.96  146.20  154.95  146.27  133.11  126.51  137.25  128.55  136.99 .131.88  143.35  132.88  139.29  137.49  135.81  132.40  161.42  151.95  160.77  149.77  0.463 KP,CP,RP  0.382,0.363 0.343  KP,AN,RP  0.416,0.363 0.353  CP,AN,RP  0.463,0.345 0.356  KP,CP,AN  0.382,0.418  RP  0.363,0.457 0.342,0.352  41  Table 4. GVRD-Ozone(n=84) Notation: Marpole(MA), Confederation Park(CP), Anmore(AN), Rocky Point Park(RP). Likelihood for dependence parameters. negative log likelihood  estimate subset  H-R  8's  r's  normal  H-R model  MA, CP  0.521  89.17  89.23  88.00  89.69  MA, A N  0.829  101.63  102.47  101.80  104.54  MA, RP  0.763  99.94  100.50  100.88  102.21  CP, A N .  0.747  98.78  99.60  95.96  101.93  CP, R P  0.680  97.23  97.72  93.19  98.95  AN, RP  0.284  61.55  61.73  72.02  64.66  MA,CP,AN  0.526,0.803  129.02  131.08  126.47  133.31  125.73  129.42  123.81  130.07  103.97  105.08  115.26  108.96  101.13  101.70  107.48  105.76  130.01  133.28  137.97  136.71  0.741 MA,CP,RP  0.526,0.738 0.676  MA,AN,RP  0.816,0.747 0.284  CP,AN,RP  0.734,0.670 0.285  MA,CP,AN  0.505,0.763  RP  0.705,0.700 0.642,0.291  42  Table 5. GVRD-Station 5 (Confederation Park). (n=64) Likelihood for dependence parameters. estimate subset  negative log likelihood  H-R  <5's  r's  normal  H - R model  1.359  127.33  126.92  125.49  126.77  so ,o  1.033  122.83  122.88  122.10  123.01  N0 ,0  0.932  120.31  120.56  118.39  120.76  1.312,1.023  181.73  180.37  176.24  179.48  S0 ,N0 2  2  2  2  3  3  S0 ,N0 ,0 2  2  3  0.930  Table 6. GVRD-Station 7 (Anmore). (n=76) Likelihood for dependence parameters. estimate subset  negative log likelihood  H-R  <5's  T'S  normal  H - R model  1.167  148.68  148.81  146.10  148.84  so ,o  1.134  148.21  148.22  145.17  148.25  N0 ,0  0.526  118.65  119.11  124.25  120.64  1.132,1.112  189.37  190.11  192.69  191.96  S0 ,N0 2  2  2  2  3  3  S0 ,N0 ,0 2  2  3  0.526  43  Table 7. GVRD-Station 9 (Rocky Point Park). (n=63) Likelihood for dependence parameters. estimate subset  negative log likelihood  H-R  tf's  r's  normal  H-R model  1.065  120.87  121.48  118.51  121.99  so ,o  1.063  120.38  121.10  118.57  122.28  N0 ,0  0.622  107.55  106.99  106.00  106.30  1.019,1.050  162.60  163.61  159.92  164.04  S0 ,N0 2  2  2  2  3  3  S0 ,N0 ,0 2  2  3  0.625  44  Chapter 4. Analysis With Missing Data In statistical analysis, it may happen by accident or some other reason that some of the observations are missing, cannot be collected, or in some other way are not obtainable. In such a situation it is obvious that the routine method is not appropriate. Now  we would like to reanalyze the data from Greater Vancouver Regional Dis-  trict with missing values among the hourly measurements. If there are missing hourly measurements, we treat the daily or weekly maxima as right-censored. For  right censored data, let x - be either the observed value or the right censored t  value, and let 1 if Xi is right-censored Si 0 if Xi is not right-censored If f(t; 0) is a parametric family of models for the iid data and if S(t; 9) is the survival function, then the likelihood function with right-censored data is L = n^J(x ,0) - S(t ,9) , 1  S i  6 i  i  i  see, for example, L'awless(1980) for details. We will use the procedure in Section 2 proposed by Joe (1990b) for fitting of a parametric multivariate distribution to independent and identically distributed p-vectors X i , . . . ,X . n  The first step is to fit generalized extreme value distribution to the uni-  variate margin separatly by maximum likelihood. That is, the j t h univariate margin is: Fj{x^ Hj,(Tj) h  =  e x p { - ( l + 7 [(x i  -  fij)/<Tj])+ }, lhi  where (z)+ = max(z,0). Assuming that the measurements are iid, the likelihood is  45  where ~F(x; 7, [i, a) = 1 - F ( x ; 7, (i, a) and f(x; 7, fi, a) = — been added for the jth  , and subscript j has  component.  After fitting the univariate margins, make the transformation  ^• = (i+7i[(^-^)/^]r  1/7>  so that original random variables become exponential random variables. The y ^ s are not treated as left censored data. Now for multivariate parameters, we take the univariate margins to be exponential.  {Y\\, Yi ),.  For bivariate data, suppose we have iid random pairs observe  (Yu, Y , Sn, 6 ),... 12  ,{Y , Y , S ,  12  nX  n2  Sij = {  ^12),and  S ),where  nl  n2  1  if Yij is left censored,  0  if Y{j is not censored,  i = 1 , . . . , n and j = 1,2. The bivariate cdf is F ( y , y ,6) x  F is the survival function and F\ and F  • • i(Yni,  2  2  = F = 1 - F - F + F, where x  2  are univariate margins of F.  2  Therefore we can get the following contribution to the likelihood of the bivariate parameter 6:  P(Y =y ,Y  2  P(Y <y ,Y  2  1  1  1  1  = y )_  dF dy dy  = y) _  dF_ _ df_ _ d~F dy dy dy  2  2  _  P(Y =y ,Y <y ) 1  1  2  2  _ P(Y <y .Y <y ) = 1  1  2  2  2  1  if  ( £ i , M =  (0,0)  2  7  2  2  8F_ _  dyi  9F _  dyi  2  = (1,0)  i2  dF, dyiif (^1,^2) =  F( ,y ,0)  if (6 S )  2  y i  if (Siu6 )  iu  =  i2  (0,1) (l,l)  In the trivariate case, the generalization is straightforward and some details are given below. and observe  Suppose we have iid random variables  (Y ,Y ,Y ,8 ,S ,S ), n  1 2  1 3  n  1 2  1 3  8^ = <  (Y\\,Yi ,Y\ ), 2  ... ,(Y i, Y , Y ),  3  n  (Y ,Y ,Y ,8 ,8 ,8 ), n l  n 2  n 3  n l  1  if Yij is left censored,  0  if Yij is uncensored, 46  n 2  n 3  where  n2  n3  i — 1,..., n and j = 1,2,3. The trivariate cdf is F(yi, y , y , 5 ) = 1 - F i ( y i ) - F ( y i ) - F ( y i ) + F ( , 2  2  2  + ^ i ( y i , y , 0) + 3  3  ^23(2/2,2/3,0)  3  12  yi  2  1  3  i2  2  - F( , y , y ,9)  where JF is survival function and F\, F , F , F , F , and F 2  y ,5)  yi  13  2 3  2  are the univariate and  bivariate margins of F. Let C i  = (Sn,Si2,Si ), then we have  2 3  3  P(Y =y ,Y  = y ,Y = y )=:  P(Y <y Y  = y ,Y  1  1  1  u  2  2  2  3  2  dF dyidy dy 3  3  2  F = y ) =dyddy3 2  3  1  P(Y  1  = ,Y Vl  < y ,Y  2  2  2  = y ,Y = y ,Y <y ) 1  2  2  3  1  2  2  3  P(Yr < y Y u  2  3  u  2  2  P{Y < ,Y 1  Vl  2  3  =  3  < y ,Y 2  3  dF dyi  d F dyidy  dF dyi dy  2  2  d F dy dy  3  dF dy\dy  2  _  2  dF* dys  _  _  dF dy  3  , dF-\3 " " dy  1  ,  1  3  dFi ?  dy  2  dF, , dF-\ 9 dyi " " dyi t  dF?3 _  1  3  dF 3 2  dy  2  2  dF_  _ a£ dy  2 3  = (0,1,0) = (0,0,1)  if C*i = (1,1,0) if  C\2  3  = (1,0,1)  2  dF-] 3 _ 9 F dyi dyi  2  3  23  dy3  < y ) = F(yi, 2/2,y ,0) 3  C\2  if C i  2  dy  +  2  if  3  2  n  2  3  P(Y =y Y <y ,Y <y ) 1  2  < y ) = dF  dy dy3  u  3  = y ,Y  2  dF dyidy  = y ) =dF  P(Y <y ,Y <y ,Y 1  =  3  = (0,0,0)  23  dyidyz 2  3  2 3  if C*i = (1,0,0)  2  2  dF = y ) = dyidy3  3  dF  d F,3  3  2  P(Y  if C i 2  if C ' i = (0,1,1) 23  if C*i  23  =  (1,1,1)  Clearly, this can be generalized to higher dimensions. Computer implementation of likelihood For models (2-7-1) to (2-7-8), the different contributions to the likelihood can easily be obtained using a symbolic manipulation software, such as Maple (Char et al., 1988, 5th edition), where derivatives can be output in Fortran format. Data analysis Here we reanalysize the Data Set 2 of Chapter 3 with the missing hourly values and censored weekly maxima. The results of the multivariate extreme value analysis are summarized in Table 8 — 10. The values of the maximum likelihood estimates are given in the left hand side of each table, and the negative log likelihoods are in the right side of each table. 47  A larger value of n means that the maxima for the station tend to be larger, a larger value of a means more spread and a larger value of 7 means that the hourly concentration has a heavier tail. Compared to the results of previous analysis, where weeks with many missing hourly measurements were treated as missing at random, we see that all parameter estimates of fi,cr and 7 in three tables are mostly slightly larger (compare Tables 3 to 5 of Joe 1990b). This is not a surprising result. After fitting the univariate margins, we obtain the transformed data which is left censored. The next step is to analyze the transformed data. The models for the multivariate exponential distribution were described before. The parameter estimates for the censored data are quite close to the estimates in Table 3 to 5 of Joe (1990b), but the relative values of the estimates are different in some cases. For trivariate case, the likelihoods for the 6's models and the r's models diverge whereas they are close in the bivariate case. A n explanation of the multivariate parts of Tables 8 to 10 is explained next for one case. For example, for (CP, SN, AN) in Table 8, the estimates of the parameters in model (2-7-3) are (1.567, 1.249), and the estimates of the parameters in (2-7-4) are (0.851, 0.507). For model (2-7-3) (respectively (2-7-4)), 1.567 (0.851) measures the dependence of the pair of stations (CP, SN), and 1.249 (0.507) measures the dependence of the pairs of (CP, A N ) and (SN, A N ) . The model (2-7-3) is better in Table 10, and the model (2-7-4) is better in Table 8 and 9. Conclusions are the same as in Joe (1990b) which means that there is no substantial difference with the missing values deleted.  48  Table 8. GVRD-Sulphur Dioxide(n=104) Notation: Confederation Park(CP), Second Narrows(SN), Anmore(AN),Rocky Point Park(RP). Likelihood for dependence parameters. subset  estimate  estimate  negative log likelihood  S's  T'S  <5's  r's  CP, SN  1.560  0.843  185.581  185.664  CP, A N  1.220  0.468  205.635  205.345  CP, R P  1.462  0.729  200.587  200.556  SN, A N  1.292  0.535  191.869  192.223  SN, R P  1.447  0.720  190.815  190.669  AN, RP  1.647  0.905  195.814  197.452  Multivariate CP,SN,AN  1.567,1.249  0.851, 0.507  284.370  283.937  CP,SN,RP  1.550,1.456  0.823, 0.674  277.888  279.263  AN,RP,CP  1.650,1.348  0.875, 0.562  291.025  293.808  AN,RP,SN  1.675,1.395  0.868, 0.596  278.128  282.235  Estimate for univariate margin parameter CP  0.313  1.857  0.834  SN  0.181  2.083  0.885  AN  0.261  3.215  2.024  RP  0.080  2.058  1.094  49  Table 9. GVRD-Nitrogen Dioxide(n=84) Notation: Kensington Park(KP), Confederation Park(CP), Anmore(AN),Rocky Point Park(RP). Likelihood for dependence parameters. estimate  estimate  <5's  T'S  8's  r's  KP, CP  2.816  2.118  161.095  161.041  KP, A N  2.673  1.965  161.941  161.944  KP, RP  2.947  2.246  154.375  154.356  CP, A N  2.168  1.482  176.638  176.104  CP, R P  3.256  2.557  144.079  143.930  AN, RP  2.968  2.272  152.358  152.171  subset  negative log likelihood  Multivariate KP,CP,AN  2.741, 2.259  2.061, 1.715  219.517  215.709  RP,KP,CP  2.945, 2.718  2.364, 2.069  194.823  189.065  RP,KP,AN  2.974, 2.597  2.226, 2.150  198.191  192.105  RP,CP,AN  2.961, 2.419  2.2.9, 1.934  204.392  199.517  Estimate for univariate margin parameter KP  0.082  5.701  2.186  CP  0.179  6.093  2.442  AN  0.086  3.750  1.760  RP  0.078  5.228  1.658  50  Table 10. GVRD-Ozone(n=84) Notation: Marpole(MA), Confederation Park(CP), Anmore(AN),Rocky Point Park(RP). Likelihood for dependence parameters. estimate  estimate  6's  r's  6'a  r's  MA, CP  2.077  1.360  143.604  143.340  MA, A N  1.726  1.012  159.720  160.279  MA, RP  1.691  0.967  159.913  160.312  CP, A N  2.174  1.467  147.824  147.753  CP, R P  2.162  1.444  144.703  144.804  AN, RP  3.616  2.917  122.432  122.622  subset  negative log likelihood  Multivariate MA,CP,AN  2.061, 1.835  1.322, 1.140  211.020  210.141  MA,CP,RP  1.948, 1.948  1.328, 1.073  207.476  210.365  AN,RP,MA  3.606, 1.759  2.868, 0.963  189.255  191.630  AN,RP,CP  3.586, 2.131  2.852, 1.439  176.798  177.085  Estimate for univariate margin parameter MA  -0.049  4.182  0.816  CP  -0.019  5.054  1.384  AN  0.047  5.250  1.595  RP  0.091  5.054  1.761  51  Chapter 5. Conclusion and Further Research This thesis has applied maximum likelihood estimation of several parametric families of multivariate extreme value distribution, including the case where some maxima are right censored. The Husler-Reiss family has more flexibility, with one parameter for each bivariate margin, than the parametric families used in Joe (1990b). However, its form is much more complicated and for diminsion p > 3, time-consuming integrations and computations of the multivariate normal cdf are needed. Also, unlike the other families, symbolic manipulation software cannot be used to obtain the density function from the survival function. The maximum likelihood estimation takes 10-20 hours on a Sun SPARCstation 330 for p =4, and is expected to take weeks for p =5 ( the programming for this case has not been completed). If both the multivariate dependence parameters and the univariate marginal parameters are simultaneously estimated, the C P U time is much more than the preceding (this also has not been done). Fortunately, it seems possible from the maximum likelihood estimates of the bivariate likelihoods to check whether the simpler models in Joe (1990b) fit the data adequately, if so, the bivariate parameter estimates must follow certain patterns. If the simpler models appear not to be adequate, then the best model, among known ones, would likely be the Hiisler-Reiss model. Further research is needed to derive families of multivariate extreme value distributions that have a dependence parameter for each bivariate margin and that have a simpler form than the Husler-Reiss model (simpler in the sense that maximum likelihood estimation becomes possible for p > 5). Some positive work in this direction has been initiated.  52  BIBLIOGRAPHY Char, Bruce W . , Geddes, K . O . , Gonnet, G . H . , Monagan, M . B . , and Watt, S.M. (1989). Maple. 5th ed., W A T C O M Publications Limited.  Coles, S. G . and Tawn, J . A . (1991). Modelling extreme multivariate events. J . R. Statist. Soc. 52, 377-392.  Galambos, J . (1987). The Asymptitic Theory of Extreme Order Statistics, second edition. Kreiger Publishing Co., Malabar, F L  Hiisler, J . and Reiss, R . - D . (1989).  Maxima of normal random vectors between  independence and complete dependence. Statistics &; Probability Letters 7 283-286.  Joe H . (1990a). Families of min-stable multivariate exponential and multivariate extrem value distributions. Statistics & Probability Letters 9 (1990) 75-81.  Joe H . (1990b). On fitting multivariate extreme value distributions, with applications to environmental data. Technical report. Department of Statistics, University of British Columbia.  Lawless J . F . (1982). Statistical Models and Methods for Lifetime data.  JOHN  W I L E Y & SONS, Inc., New York.  Resnick, S.I. (1987). Extreme Values, Regular Variation and Point Processes. New York: Springer-verlag.  53  Schervish, J . (1984).  Multivariate normal probabilities with error bound.  Appl.  Stat.33 81- 94.  Tawn J . A . (1988). Bivariate extreme value theory, models and estimation, Biometrika 75, 397-415.  Tawn J . A . (1990). Modelling multivariate extreme value distributions. Biometrika 77, 245-253.  54  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Country Views Downloads
China 13 0
United States 3 0
Germany 1 0
United Kingdom 1 0
Canada 1 0
France 1 0
City Views Downloads
Shenzhen 9 0
Beijing 3 0
Ashburn 3 0
Unknown 2 11
Toronto 1 0
Chemnitz 1 0
Hangzhou 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}
Download Stats

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0302381/manifest

Comment

Related Items