Rebuttal of “Problems in the extreme value analysis”

Rebuttal of “Problems in the extreme value analysis”

Structural Safety 34 (2012) 418–423 Contents lists available at SciVerse ScienceDirect Structural Safety journal homepage: www.elsevier.com/locate/s...

507KB Sizes 2 Downloads 64 Views

Structural Safety 34 (2012) 418–423

Contents lists available at SciVerse ScienceDirect

Structural Safety journal homepage: www.elsevier.com/locate/strusafe

Rebuttal of ‘‘Problems in the extreme value analysis’’ Nicholas J. Cook ⇑ Department of Civil Engineering, University of Birmingham, UK

a r t i c l e

i n f o

Article history: Received 8 October 2010 Received in revised form 28 June 2011 Accepted 11 August 2011 Available online 15 September 2011 Keywords: Extremes Extreme value analysis Probabilistic design Structural safety Return period Plotting positions Risk analysis

a b s t r a c t This paper rebuts claims published in this journal [33] as well as elsewhere [31,32,34]: (a) that improved methodologies for extreme value analysis (EVA) developed over the past 60 years are invalid; (b) that EVA methodologies should revert to the status quo ante 1939; and (c) that, consequently, all regulations and codes of practice for extreme winds should be reassessed. This paper rebuts these claims and shows current EVA methodologies to be valid. The paper also shows that uncertainty due to sampling error, viz. how well a single observed sample represents the random process sampled, dominates over the choice of methodology. Ó 2011 Elsevier Ltd. All rights reserved.

1. Introduction The claim is made in Makkonen [33] – hereafter ‘‘M’’, for brevity – and elsewhere [31,32] that Weibull’s [41] estimator is uniquely the only valid estimator for the probability of ranked samples, thus invalidating all improvements in extreme value analysis methods developed since 1939, so that the many estimates of weather-related risks need to be re-evaluated and the building codes and other related regulations updated. Seeking to overturn 60 years of theoretical and practical development from two generations of eminent statisticians is a bold venture, and one would think it would have been accompanied by a rigorous and comprehensive proof. But the claim in ‘‘M’’ rests solely on untested assertions that the Weibull estimator is ‘‘unique’’ and is ‘‘exact’’. A review by Cook [9] of the first paper of the series, [31], used rigorous statistical proofs to show this is not true. This paper rebuts the further claims made in ‘‘M’’, demonstrating that the unbiased methods developed since 1939 are valid and showing that the best of these methods provide a significant improvement in accuracy. It also shows that, in practice, the uncertainty due to sampling error dominates over the choice of plottingposition methodology.

distribution of sampling error. A consistent notation is required to distinguish between the characteristics of a random process and the characteristics of samples drawn from this process. Consider a random process of the variate x that exists infinitely into the past and into the future, where the corresponding probability, P, has a cumulative distribution function (CDF), P(x). The CDF is a single-valued function which is reversible, i.e. a mapping function. Hence:

P ¼ PðxÞ and x ¼ P1 ðPÞ

Analysis of order statistics uses independent samples ranked in ascending, or descending order of value. Here the subscript m:N is used to indicate the mth out of N values of the variate x, ranked in ascending order of value, m = 1, . . . , N, i.e. the mth lowest value. As Pm:N = P(xm:N) is a single valued function of xm:N, with values Pm:N lying between 0 and 1, any general function of the variate, G(xm:N), can be substituted by G(Pm:N) and integration limits for dx can be replaced by 0 and 1 for dP. The principal issue in this paper involves the interpretation of ensemble averages, so the notation h in is used to indicate the ensemble average of n independent samples:

2. Establishing a consistent notation Cook [9] showed that the principal reason for the erroneous claim is a failure to recognise that every estimate has an associated ⇑ Tel.: +44 1727 764719; fax: +44 1727 764720. E-mail addresses: [email protected], [email protected] 0167-4730/$ - see front matter Ó 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.strusafe.2011.08.002

ð1Þ

hxin ¼

n 1X xi n i¼1

ð2Þ

Thus hxm:Ni1 indicates the ensemble average of xm:N for all possible samples, n = 1, which is conventionally called the expectation, E(xm:N).

N.J. Cook / Structural Safety 34 (2012) 418–423

419

3. Expectation and standard error

5. Reduced variate

The usual aim of extreme value analysis (EVA) is to estimate the properties of the random process, especially of P(x), to some defined degree of accuracy, using observations of x made over a finite sampling period to determine the corresponding observed sample frequency F(x). If it were possible to acquire an infinite number of independent samples of size N, then:

The claim in ‘‘M’’ that the Weibull estimator is ‘‘unique’’ rests on the presumption that the one and only purpose of EVA is to determine the probability of a given sample value. But in most cases, the principal purpose is to obtain the best unbiased estimate of the variate for a given design value of probability. In EVA, this usually means fitting the expected EV model to the observations so that the model may be extrapolated past the observations into the upper tail. The Gumbel reduced variate g is defined as:

hFðxm:N Þi1  EðPm:N Þ  Pðxm:N Þ

ð3Þ

because the ensemble mean frequency of an infinite number of independent samples is, by definition, identical to the CDF of the random process. In practice, EVA is performed on a single observed sample, and the difference between this single sample and the expectation is taken as the sampling error. For ranked independent observations of a random variate, the theoretical expectation for any function G(Pm:N) is defined by the Binomial:

EðGðPm:N ÞÞ ¼

N! ðm  1Þ!ðN  mÞ!

Z

1

GðPm:N ÞPm1 ½1  PNm dP

ð4Þ

0

It is important to note here that (4) is universally applicable to all functions of Pm:N, including probability itself, i.e. including G(Pm:N) = Pm:N. The corresponding standard deviation is given by:

rðGðPm:N ÞÞ ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Eð½GðPm:N Þ2 Þ  ½EðGðPm:N ÞÞ2

ð5Þ

For an ensemble average of n samples, the standard error en is given by:

en ðGðPm:N ÞÞ ¼

rðGðPm:N ÞÞ pffiffiffi n

ð6Þ

Hence, the standard error is equal to the standard deviation for the usual case of a single sample. 4. Weibull’s estimator The Weibull [41] estimator for Pm:N is founded on the generally accepted principle that the best unbiased estimate of any probabilistic variate is its expectation. E(Pm:N) is obtained by substituting G(Pm:N) = Pm:N into (4) [14], Eq. 2.1.3(1); [5], Eq. 2.41) and evaluating, leading to:

EðP m:N Þ ¼ hFðxm:N Þi1 ¼

m Nþ1

ð7Þ

which corresponds to the mean frequency of the mth ranked value. The standard deviation of Pm:N evaluates through (2) and (3) as:

rðPm:N Þ ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi mðN  m þ 1Þ EðP 2m:N Þ  ½EðPm:N Þ2 ¼ ðN þ 1Þ2 ðN þ 2Þ

m  Nþ1

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi mðN  m þ 1Þ ðN þ 1Þ2 ðN þ 2Þ

ð10Þ

This was originally derived for application with the asymptotic FT1 distribution:

P ¼ expð expðgÞÞ

ð11Þ

but may be viewed simply as the standard normalisation for an EV variate, in the same way that (x  l)/r is the standard normalisation for a Normal variable. In any case, as g and x are linearly related, the best unbiased estimate of x is given by the expectation of g(x). For an asymptotic FT1 distribution, the expectation of g is given by putting G(Pm:N) = ln(ln(1  P)) into (4). Unfortunately, this cannot be directly evaluated, requiring a numerical solution. In practice, there is a different solution for every value of m and N [17], but the estimator of Gringorten [13] is a very good general approximation:

PðEðgm:N ÞÞ 

m  0:44 N þ 0:12

ð12Þ

(Note: Although McClung and Mears’ [35] estimator, P = (m  0.4)/ N, appears very similar, it is a less accurate approximation – see Table 1, later.) Unfortunately, Gringorten does not derive the standard error associated with (12). Strictly, E(gm:N) relates only to the asymptotic FT1 distribution but, in practice, it is applicable to any distribution that plots as a gentle curve on Gumbel1 axes [21]. 6. Rebuttal of the claims in ‘‘M’’ As the series of papers culminating in ‘‘M’’, [31–33], consist largely of a series of untested assertions, one option is to test each assertion in turn, as in Cook [9] for [31]. This paper concentrates only rebutting just those issues in ‘‘M’’ and in the response to the earlier rebuttal [34] that are most significant to structural safety, viz: the assertion that the Weibull estimator is ‘‘unique’’ and ‘‘exact’’ and the consequences for design codes and regulations stemming from this. 6.1. The problem with the plotting position

ð8Þ

([14], Eq. 2.1.3(6); [5], Eq. 2.42). Weibull’s estimator is therefore revealed as being one particular application of the universal expression for expectation E( ), namely the solution applicable to the mean frequency of Pm:N. The standard error in applying this to the frequency of a single sample is r(Pm:N). Hence Eq. (12) for Pm:N in ‘‘M’’ is not exact, and should be expanded to:

Pm:N ¼

g ¼ aðx  UÞ

ð9Þ

The principal issue in Section 3 of ‘‘M’’ is the assertion that ‘‘there exists only one correct plotting formula’’, viz. that given by (7). In practice, no estimator for a sample of data is ever perfect and the ‘‘best’’ estimator in any context depends on variable being estimated the purpose to which this estimate is put. Plotting positions that depend on the expected model distribution are conventionally used to test whether that model is appropriate. On axes that linearise the model, a linear fit is therefore indicative of the good fit to the model. There is no ‘‘circular deduction’’ here, 1 Although the practice of plotting the variable as ordinate and ln(ln(P)) as abscissa is traditionally called a ‘‘Gumbel plot’’, or ‘‘Gumbel axes’’, the first use of these axes is attributed by Gumbel [14, P 176] to Powell [37].

420

N.J. Cook / Structural Safety 34 (2012) 418–423 Table 1 The design reduced variate g(P = 0.98) for an un-weighted least mean square fit to N = 21 extremes, for the methods in Makkonen [33] (re-sorted chronologically) and the percentage bias compared with the datum of E(g). Estimator for Pm:N

Reference

g(P = 0.98)

Bias in g (%)

(m  0.5)/N m/N m/(N + 1) (m  0.31)/(N + 0.38) (m  0.375)/(N + 0.25) (3m  1)/(3N + 1) (m  0.44)/(N + 0.12) exp(exp((0.577 + 3ln2(2m/(N + 1)  1)))) (m  0.35)/N (m  0.4)/N Numerical method

Hazen [23] Anon [1] Weibull [41] Beard [3] Blom [4] Tukey [40] Gringorten [13] Barnett [2] Landwehr et al. [28] McClung and Mears [35] Harris [16]

3.822 3.578 4.303 4.024 3.959 4.001 3.890 3.879 3.642 3.708 3.902

2.06 8.30 +10.28 +3.14 +1.46 +2.54 0.32 0.59 6.65 4.97 0

just the standard scientific method: suggest a model ? test the model ? accept/reject the model. Table 1 and Fig. 1 in ‘‘M’’ were built by adopting the Weibull estimator for E(P) as the datum estimator, but the purpose of the table is to compare the predicted return period, R, of the highest rank m = N = 21. The correct datum is therefore the probability corresponding to the expectation of R, which is the solution of (4) for G(P) = 1/(1  P), viz. P(E(Rm:N)) = m/N, leading to RN:N = 1 as the expectation. For structural safety purposes, the plotting position which gives the best unbiased estimate of the variate for a given design probability is preferred. Gringorten’s [13] estimator and Harris’ [17] numerical method give almost identical values of this. It is not a problem that these produce ‘‘a parameter that cannot be re-transformed to the non-exceedence probability of the variable’’2 because this non-concomitance is expected, indeed is required to eliminate bias in the value of the variate. As long as N is finite, i.e. for all practical analyses of real data, there is never one single solution, but there will be a best unbiased solution for each purpose which will not generally be concomitant with other solutions: E(P) – P(E(R)) – P(E(g)). Similarly Fig. 1 in ‘‘M’’ also assumes the Weibull estimator as the datum estimator. The usual purpose of a Gumbel plot is to obtain the fitted parameters of a FT1 distribution to the values, i.e. to give the unbiased estimate, E(x(P)), of the variate, x, for any given

2

Here Makkonen means to say ‘‘cannot be re-transformed after ensemble averaging’’ as it is the averaging process that destroys concomitance.

η

Fig. 1. Makkonen’s Fig. 1 data for the top six ranked wind speeds of sample size N = 21, re-plotted on Gumbel axes for Weibull, Gringorten and return period estimators. See key on figure to identify curves.

risk. The data from Fig. 1 in ‘‘M’’, which represent the top 6 ranks out of 21 sampled values, have been re-plotted in Fig. 1 here, using the standard Gumbel axes to represent the asymptotic FT1 model as a straight line. The corresponding distributions of sampling error for these top 6 ranks, obtained to a precision of 1% by using 104 Bootstrapped samples, are shown in Fig. 2. The Weibull plotting position, marked for each curve in Fig. 2 by the lines tagged ‘‘W’’, is seen to correspond with the mode of the sampling error. The Weibull estimator would therefore be appropriate when the most likely value is required. The Gringorten plotting position, marked for each curve in Fig. 2 by the lines tagged ‘‘G’’, is seen to represent the mean of the sampling error. It is appropriate when fitting a model to the observations using methodologies that require the mean value (e.g. least mean square or maximum likelihood methods) to obtain an unbiased fit. The curve appropriate for the return period datum, E(R), is included in Fig. 1 to demonstrate that the highest rank cannot be plotted because E(RN:N) = 1. The reciprocal nature of return period makes it a very poor datum for comparing errors in the variate, because equal increments of error produce increasingly greater increments in return period. In any case, as the largest observation is unusable, there is little point using return period for comparison and so Table 1 in ‘‘M’’ has little practical relevance. A better basis of comparison is the reduced variate for the typical design risk, g(P = 0.98), because this is linearly related to the design value. Table 1 in this paper compares the design reduced variate for an unweighted least mean square fit to N = 21 extremes, for each of the methods in Table 1 of Makkonen [33] (re-sorted chronologically),

Fig. 2. Density function from 10,000 Bootstrap samples of the reduced variate for the top 6 ranks of sample size N = 21. See key to identify ranks. ‘‘W’’ and ‘‘G’’ mark the Weibull [41] and Gringorten [13] plotting positions, respectively.

N.J. Cook / Structural Safety 34 (2012) 418–423

together with the bias error for the datum of E(g). The largest overestimate (+10.28%) is given by the Weibull [41] estimator, while the largest underestimate (6.65%) is given by the Landwehr et al. [28] estimator. For the two methods that adopt E(g) as the datum the error is theoretically zero but, in practice, the Gringorten [13] estimator is a general approximation for all N and exhibits as small error, whereas the Harris [16] method is exact to the computed accuracy. The percentage errors in the variate, x, are considerably smaller than the corresponding percentage errors for the reduced variate, g, because of the diluting effect of the annual modal value. For example, for annual maximum wind speeds, V, in the UK the +10.28% bias for g(P = 0.98) translates to a +2.89% bias for V(P = 0.98). Note that the mode and dispersion values used in Fig. 1 in ‘‘M’’ are far from typical of temperate zone wind speeds: aU (V) = 10.46  0.175 = 1.8 in ‘‘M’’, whereas typically aU (V)  10. Also, only the least reliable top 6 ranks are fitted in ‘‘M’’ instead of the full 21. These two actions grossly exaggerate the errors in wind speed and dynamic pressure quoted in the caption to this figure. Fig. 2, here, demonstrates that the uncertainty due to sampling error dominates over the choice of plotting position method. 6.2. The solution to the plotting positions In Section 4 of ‘‘M’’ the route to obtaining the Weibull estimator for Pm:N is correctly derived and demonstrated to represent the ensemble mean for an infinite number of samples. But the standard error that results from applying this result to a single sample is not acknowledged. Indeed, in response to Cook [9], Makkonen [34] reiterates the claim that the result is ‘‘exact’’ because ‘‘there is no sample error associated’’ – confusing the derivation of the Weibull estimator from all possible samples (which is exact) with its application to a single sample (which incurs the corresponding standard error given by (8)). It is important to recognise that there are two limiting processes involved: 1. Expectation, E( ), and standard deviation, r( ), are derived in the limit as the number of samples n ? 1. In practical observations n = 1. 2. E(P) – P(E(R)) – P(E(g)), have different dependence on the sample size N, and converge together as N ? 1. However, samples sizes are often quite small in practice, so a significant difference is to be expected unless the variables are linearly related. For the assertion that ‘‘all suggested plotting formulas other than the Weibull formula are incorrect’’ to be true, the principle of expectation must apply exclusively to Pm:N, instead of being universal. This begs the question of when other valid plotting formulas are appropriate. The arguments in this section all relate to the exclusive assessment of E(P), with no admission that other parameters are more appropriate for specific purposes, i.e. E(g) to give the expectation of the variate for engineering design. 6.3. Errors in determining the plotting positions Contrary to the claim in 5.2 in ‘‘M’’, the choice of the ordinate scale on the standard Gumbel axes is not ‘‘arbitrary’’, but is chosen so that the asymptotic FT1 model appears as a straight line. In this context a good linear fit is far from being an ‘‘illusion’’ since a straight line corresponds to a good fit to the model. The estimator corresponding to E(gm:N) is required when mean fitting methods are used if the fit is to be unbiased. The comment by Castillo [5] that ‘‘the plotting position formulas can affect the linear trend of the cumulative probability distribution so that a careful selection must be made’’ does not illustrate ‘‘confusion’’ on his part, but rather that he recognises that there is a choice which depends on the purpose of the analysis and it is important

421

to select the correct one. If there is ‘‘confusion’’, it is certainly not in the minds of Lieblein [29], Gringorten [13], Castillo [5], Harris [17,19], Cook [6], Cook and Harris [8] and the other authors of the unbiased methods that are dismissed as being ‘‘wrong’’. Each of these authors clearly states the aims of their methodology and executes it rigorously. 6.4. The applicability of the asymptotic extreme value theory Like many other authors, Makkonen attributes the success of conventional EVA of wind speed data to the fact that there is the large number of 8766 hourly values in each year. It has been known for a very long time [11] that there are only around 150 independent wind events in each year. Simiu and Heckert [39] show that 2-day separated wind maxima are not quite independent, while 4-day separated wind maxima probably are. Harris [20] recently showed that the correlation time for hourly-mean wind speeds is about 22 h: and the usual assumption of three times the correlation time also gives a good match. It is only recently that methods which use the general or ‘‘penultimate’’ FT1 model to avoid the problems of asymptotic convergence have been advocated [8] [21]. The final paragraph of Section 6.1 in ‘‘M’’, implies that un-converged distributions are typically different in the body than in the upper tail and typically exhibit an abrupt bend in the curve. Neither is true: effects of non-convergence decrease smoothly through the range and there are no abrupt changes. A ‘‘bend’’ in the data is indicative of two or more physical mechanisms and analysis of such data should be performed using a multi-mechanism model, as in the methodologies of Gomes and Vickery [12] or Cook et al. [7]. 6.5. The fitting method Section 7 in ‘‘M’’ reviews the subject of fitting weights. The common approach of inversely weighting by the variance (for method of least squares), or standard deviation (for maximum likelihood) of the sampling error, does not ‘‘effectively neglect’’ the upper tail, but weights each point according to its statistical reliability. In addressing the dominant mechanism in a mixed climate by fitting only to the upper tail of the observations, the weights for those observations not used should be discarded and the weights used should be re-normalised to unity: then each observation will still be appropriately weighted with respect to the others. A much better approach for mixed climates is to separate out the observations by causal mechanism as demonstrated by Gomes and Vickery [12], Cook et al. [7] and, more recently, Lombardo et al. [30]. 6.6. Discussion and conclusions In Section 8 of ‘‘M’’, the choice of ordinate, abscissa and the error to minimise is inconsistent. Castillo [5] and Harris [17] favour plotting with transposed Gumbel axes: ln(ln(Pm:N)) as abscissa and xm:N as ordinate, precisely because Pm:N is a random variable and sampling error dominates over the measurement error3 associated with xm:N. The natural fit on these axes is to minimise errors in ln(ln(Pm:N)). This author favours the standard Gumbel axes because they are familiar to most practitioners, and mode & dispersion correspond to intercept & slope, respectively: but also favours minimising errors in ln(ln(Pm:N)). When fitting to a linear model the choice is irrelevant because xm:N and ln(ln(Pm:N)) remain concomitant. 3 For most of the archived records of wind speed, the measurement resolution of anemometers has been 1 kt (0.5 m/s) which is why records are kept in integer values of knots (WMO standard).

422

N.J. Cook / Structural Safety 34 (2012) 418–423

The final comment in Section 8 of ‘‘M’’ is that [17] method ‘‘is still being promoted in scientific literature’’ – but this is because it is rigorously correct and is among the best current methods for obtaining unbiased estimates of the variate.

7. The ‘‘unique & exact’’ fallacy The central argument in ‘‘M’’ used to condemn all developments in EVA methodology over the past 60 years and on which the various assertions and prohibitions are founded is the fallacy that the Weibull estimator [41] is ‘‘unique’’ and ‘‘exact’’. The response [34] to the rebuttal [9] of his first paper [31] restates these assertions: ‘‘the Weibull formula provides the exact distribution-free probability P of non-exceedance by a random observation drawn from the population. There is no sampling error related to P. Hence, the use of ‘‘estimators’’ for plotting positions is both unnecessary and misleading.’’ Also: ‘‘When considering order-ranked data, the sampling error is related solely to the observed variable x, not the probability of nonexceedance’’. Sampling error relates to how well the single observed sample represents the random process that is sampled – the Weibull estimator is derived for P, i.e. for n = 1, but is applied to F, i.e. to n = 1. Castillo [5] refers to this under the heading ‘‘Order statistics of one sample from a uniform parent’’. Error in wind speed values comes solely from the accuracy of measurement and precision of archiving, and this is around ±0.5 kt (±0.25 m/s), as noted earlier. The claim that sampling error relates to the value of the observed variable, and not the plotting position, is easily demonstrated to be false. Fig. 3 shows the Gumbel plot, using the Gringorten [13] plotting positions, of annual maximum gust speeds for three observation periods, each of which contains the same largest observed value. The position of the maximum gust speed (88 kt) is labelled ‘‘A’’ for the 10-year epoch, ‘‘B’’ for the 30-year epoch and ‘‘C’’ for the 70-year epoch. Increasing the epoch length results in a reduction in sampling error and a revised estimate of the plotting position of this observation but, obviously, makes no change to its value. Also shown in Fig. 3 are the 5– 95% confidence limits obtained from 104 Bootstrap samples for N = 10, 30 and 70.

Fig. 3. Gumbel plot of observed annual maximum gust speeds using the Gringorten [13] plotting positions. The open square symbols are for the 10-year epoch, 1960– 1969 and the highest observation of 88 kt is labelled ‘‘A’’. The filled circles are observations for the 30-year epoch, 1940–1969, and the plotting position for the 88 kt observation moves to ‘‘B’’. The open circles are observations for the 70-year epoch, 1940–2009, and the plotting position for the 88 kt observation moves to ‘‘C’’. This demonstrates that sampling error manifests in uncertainty of the plotting position and is not due ‘‘solely to the observed variable’’ [33].

8. Summary and conclusions The prescriptive approach in ‘‘M’’ is based on the fundamental fallacy that the Weibull estimator is ‘‘unique’’ and ‘‘exact’’, i.e. there is no valid alternative and the result exhibits no sampling error. This fundamental fallacy is used to condemn all developments in EVA methodology over the past 60 years, using quite strong terms – ‘‘error/erroneous’’ (20 times), ‘‘improper’’ (7), ‘‘misleading’’ (5), ‘‘wrong’’ (4), ‘‘misuse’’ (2) and ‘‘illogical’’ (1). The ‘‘hit list’’ of authors whose work is ‘‘invalidated’’ by this fallacy is quite impressive: Barnett [2], Blom [4], Castillo [5], Cunnane [10], Gringorten [13], Harter [22], Haktanir [15], Harris [16–18], Kimball [27], Hosking et al. [24], Jones [25], Jordaan [26], Landwehr et al. [28], McClung and Mears [35], Naess [36], Rasmussen and Gautam [38], Whalen et al. [42,43] and Yu and Huang [44]. The impact that choice of plotting position methodology makes on the result is greatly exaggerated in ‘‘M’’. The choice only significantly affects the top few ranks, where the sampling error is greatest. Bias is minimised by weighting the fit inversely to the variance of the plotting position, giving more weight to the body than to the tails of the distribution. However the bias produced by the use of the Weibull estimator may make the top few ranks appear to arise from a second mechanism, when no such second mechanism exists. This paper and a previous paper [9] use rigorous theory and demonstrations to demonstrate that the prescriptive approach to EVA advocated in ‘‘M’’ is based on a fundamental fallacy. This fallacy is used to justify a series of assertions that are diametrically opposed to current best knowledge and practice. Fortunately, these assertions have been easily disproved. The improvements to extreme value analysis methods made over the past 60 years are valid; it is not necessary or desirable to revert to the status quo ante 1939; and the current generation of building codes and regulations do not require urgent updating. Acknowledgement Thanks are due to Mr. R.I. Harris for many fruitful discussions and comments on the manuscript. References [1] Anon. Flow in California streams. Calif Dep Public Works, Div Eng Irrig, Bulletin 5; 1923. [2] Barnett V. Probability plotting methods and order statistics. Appl Statist 1975;24:95–108. [3] Beard LR. Statistical analysis in hydrology. Trans Am Soc Civ Eng 1943;3108:1110–60. [4] Blom G. Statistical estimates and transformed beta-variables. John Wiley; 1958. [5] Castillo E. Extreme value theory in engineering. Academic Press; 1988 [389pp.]. [6] Cook NJ. Towards better estimation of extreme winds. J Wind Eng Ind Aerodyn 1982;9:295–323. [7] Cook NJ, Harris RI, Whiting RJ. Extreme wind speeds in mixed climates revisited. J Wind Eng Ind Aerodyn 2003;91:403–22. [8] Cook NJ, Harris RI. Exact and general FT1 penultimate distributions of extreme wind speeds drawn from a tail-equivalent Weibull parents. Struct Safety 2004;26:391–420. [9] Cook NJ. Comments on ‘‘Plotting positions in extreme value analysis’’ (The role of sampling error in extreme value analysis). J Appl Meteorol Climatol 2011;50(1):255–66. [10] Cunnane C. Unbiased plotting positions – a review. J Hydrol 1978;37:205–22. [11] Davenport AG. Note on the distribution of the largest value of a random fluctuation with application to gust loading. Proc Inst Civil Eng 1964;28:187–96. [12] Gomes L, Vickery BJ. Extreme wind speeds in mixed climates. J Wind Eng Ind Aerodyn 1977;2:331–44. [13] Gringorten II. A plotting rule for extreme probability paper. J Geophys Res 1963;68:813–4. [14] Gumbel EJ. Statistics of extremes. Columbia Univ. Press; 1958. [15] Haktanir T. Self-determined probability-weighed moments method and its application to various distributions. J Hydrol 1997;194:180–200. [16] Harris RI. Gumbel re-visited – a new look at extreme value statistics applied to wind speeds. J Wind Eng Ind Aerodyn 1996;59:1–22.

N.J. Cook / Structural Safety 34 (2012) 418–423 [17] Harris RI. Improvements to the method of independent storms. J Wind Eng Ind Aerodyn 1999;80:1–30. [18] Harris RI. Control curves for extreme value methods. J Wind Eng Ind Aerodyn 2000;88:119–31. [19] Harris RI. Extreme value analysis of epoch maxima, convergence and choice of asymptote. J Wind Eng Ind Aerodyn 2004;92:897–918. [20] Harris RI. The macrometeorological spectrum – a preliminary study. J Wind Eng Ind Aerodyn 2008;96:2294–307. [21] Harris RI. XIMIS, a penultimate extreme value method suitable for all types of wind climate. J Wind Eng Ind Aerodyn 2009;97:271–86. [22] Harter HL. Another look at plotting positions. Comm Stat – Theor Meth 1984;13:1613–33. [23] Hazen A. Storage to be provided in impounding reservoirs for municipal water supply. Trans Am Soc Civ Eng 1914;77:1547–50. [24] Hosking JR, Wallis JR, Wood EF. Estimation of the generalized extreme-value distribution by the method of probability weighted moments. Technometrics 1985;27:251–61. [25] Jones DA. Plotting positions via maximum-likelihood for a non-standard situation. Hydrol Earth Syst Sci 1997;1:357–66. [26] Jordaan I. Decisions under uncertainty. Cambridge Univ Press; 2005. [27] Kimball BF. On the choice of plotting positions on probability paper. J Am Stat Assoc 1960;55:546–60. [28] Landwehr JM, Matalas NC, Wallis JR. Probability weighed moments compared with some traditional techniques in estimating Gumbel parameters and quantiles. Water Resour Res 1979;15:1055–64. [29] Lieblein J. Efficient methods of extreme value methodology. Report NBSIR 74602. Washington: National Bureau of Standards; 1974. [30] Lombardo FT, Main JA, Simiu E. Automated extraction and classification of thunderstorm and non-thunderstorm wind data for extreme-value analysis. J Wind Eng Ind Aerodyn 2009;97:120–31.

423

[31] Makkonen L. Plotting positions in extreme value analysis. J Appl Meteorol Climatol 2006;45:334–40. [32] Makkonen L. Bringing closure to the plotting position controversy. Comm Stat – Theor Meth 2008;37(3):460–7. [33] Makkonen L. Problems in the extreme value analysis. Struct Safety 2008;30:405–19. [34] Makkonen L. Reply to ‘‘Comments on ‘‘Plotting positions in extreme value analysis’’ by N. Cook. J Appl Meteorol Climatol 2011;50:267–70. [35] McClung DM, Mears AI. Extreme value prediction of snow avalanche runout. Cold Reg Sci Technol 1991;19:163–75. [36] Naess A. Estimation of long return period design values for wind speeds. J Eng Mech 1998;124:252–9. [37] Powell RW. A simple method of estimating flood frequencies. Civ Eng 1942;14:105. [38] Rasmussen PF, Gautam N. Alternative PWM-estimators of the Gumbel distribution. J Hydrol 2003;280:265–71. [39] Simiu E, Heckert NA. Extreme wind distribution tails: a ‘peaks over threshold’ approach. J Struct Eng 1996;122(5):539–47. [40] Tukey JW. The future of data analysis. Ann Math Stat 1962;33:1–67. [41] Weibull W. A statistical theory of strength of materials. Ing Vetenskaps Ak Handl 1939;151:1–47. [42] Whalen TM, Savage GT, Jeong GD. The method of self-determined probability weighted moments revisited. J Hydrol 2002;268:177–91. [43] Whalen TM, Savage GT, Jeong GD. An evaluation of the self-determined probability-weighted moment method for estimating extreme wind speeds. J Wind Eng Ind Aerodyn 2004;92:219–39. [44] Yu G-H, Huang C-C. A distribution free plotting position. Stoch Environ Res Risk Assess 2001;15:462–76.