JOURNAL
OF GEOPHYSICAL
RESEARCH,
VOL. 103, NO. D14, PAGES 17,149-17,161, JULY 27, 1998
Factors affecting the detection of trends: Statistical considerations and applications to environmental data Elizabeth C. Weatherhead,• Gregory C. Reinsel,2 George C. Tiao,3 Xiao-Li Meng,4 DongseokChoi,4 Wai-Kwong Cheang,2 Teddie Keller? John DeLuisi, 6 Donald J. Wuebbles,? James B. Kerr, 8 Alvin J. Miller, 9 Samuel J. Oltmans, •ø and John E. Frederick, TM Abstract. Detection of long-term, linear trends is affected by a number of factors, includingthe size of trend to be detected,the time spanof availabledata, and the magnitudeof variability and autocorrelationof the noise in the data. The number of years of data necessaryto detect a trend is stronglydependenton, and increaseswith, the
magnitude of variance (o-2•)andautocorrelation coefficient (qb)of thenoise.Fora typical rangeof valuesof o-2•and4>thenumberof yearsof dataneededto detecta trendof 5%/decadecan vary from -10 to >20 years,implyingthat in choosingsitesto detect trends somelocationsare likely to be more efficient and cost-effectivethan others. Additionally, someenvironmentalvariablesallow for an earlier detection of trends than other variablesbecauseof their low variability and autocorrelation.The detectionof trends can be confoundedwhen suddenchangesoccur in the data, suchas when an instrumentis changedor a volcano erupts. Suddenlevel shiftsin data sets,whether due to artificial sources,suchas changesin instrumentationor site location,or natural sources, suchas volcaniceruptionsor local changesto the environment,can stronglyimpact the number of years necessaryto detect a given trend, increasingthe number of years by as much as 50% or more. This paper providesformulae for estimatingthe number of years necessaryto detect trends,alongwith the estimatesof the impact of interventionson trend detection.The uncertaintyassociatedwith these estimatesis also explored.The resultspresentedare relevantfor a variety of practicaldecisionsin managinga monitoring station,suchas whether to move an instrument,changemonitoringprotocolsin the middle of a long-termmonitoringprogram, or try to reduce uncertaintyin the measurements by improvedcalibrationtechniques.The resultsare also usefulfor establishingreasonableexpectationsfor trend detectionand can be helpful in selecting sitesand environmentalvariablesfor the detectionof trends.An important implicationof theseresultsis that it will take severaldecadesof high-qualitydata to detect the trends likely to occur in nature. have attempted to detect long-term trendsin geophysicalvariablessuchas atmosphericozone [Reinselet al., 1994;Stolarski The impact of human interventionin a changingenviron- et al., 1991, 1992], stratospherictemperature [Miller et al., ment hasbroughtabout increasedconcernfor detectingtrends 1992], and ultraviolet (UV) radiation at the Earth's surface in various types of environmentaldata. A variety of studies [Scottoet al., 1988; Weatherheadet al., 1997]. These studies have revealed that detectionof trendsis difficult and requires 1Cooperative Institutefor Research in theEnvironmental Sciences, statisticaltechniquesthat take into accountsome of the realUniversity of Colorado, Boulder. isticproblemswhich frequentlyoccurwith the measurementof 2Department of Statistics, University of Wisconsin, Madison. geophysicaldata. As political decisionsand future scientific 3Graduate Schoolof Business, University of Chicago, Chicago, Illi- efforts may be based on the results of environmental trend nois. 4Department of Statistics, University of Chicago, Chicago, Illinois. studies,achievingthe most accurate trend estimatesin the SNational Centerfor Atmospheric Research, Boulder,Colorado. shortesttime period is important to ensuringthat appropriate 6NOAAAir Resources Laboratory, Boulder,Colorado. actions are taken. Factors affecting the detection of linear 7Department ofAtmospheric Sciences, University of Illinois,Cham- trends are outlined in this paper, and examplesof trend depaign. 8Atmospheric Environment Service, Downsview, Ontario,Canada. tection analysesof data are highlighted.The resultsshowthat for most expectedenvironmentalchanges,severaldecadesof 9NOAANationalWeatherService, Washington, D.C. IøNOAAClimateMonitoring andDiagnostics Laboratory, Boulder, high-qualitydata will be neededbefore suchchangeswill be Colorado. detectable.The resultsalsoshowthat certain decisionsregard11Department of Geophysical Sciences, University of Chicago, Chi- ing site location, instrumentmaintenance,and calibration can cago,Illinois. influencethe accuracyof trend estimationand hencedetection 1.
Introduction
Copyright1998 by the American GeophysicalUnion.
of trends.
Paper number 98JD00995.
Statisticalcriteria for detectinglinear trendswere presented by Tiao et al. [1990], where it was shownthat the precisionof
0148-0227/98/98JD-00995509.00
17,149
17,150
WEATHERHEAD
ET AL.: TREND DETECTION
trend estimatesis strongly influenced by the variability and autocorrelationof the underlyingnoiseprocess.The precision of a trend estimatein turn directly determinesthe number of yearsof data requiredto detecta trend of givenmagnitudeor the magnitude of trend that can be detectedwith a given number of years of data. This paper exploresmore fully, in both a theoreticaland an applied sense,how variousfactors affect the ability to detect a trend through their influenceon the precisionof the trend estimate.The uncertaintyof estimatesfor trend detectabilityis alsoexplored.Section2 of this paperaddresses the abilityto detecta trend in a singledataset and how the number of yearsof data neededto detect a given trend is dependenton the magnitudeof and autocorrelation in the naturalvariation of the data, aswell as how theseparameters can vary significantlyamong environmentalparameters as well as from location
to location.
Section 3 addresses the
temperature on one day is often associatedwith higher than normal temperatureon the next day. This positiveautocorrelation
tends
to confound
with
a linear
trend
and therefore
increasesthe length of time required to detect a given trend. 2.1.
Basic Statistical Modeling
In manystatisticaltrend analysesof monthlytime seriesdata on a geophysicalvariable, such as total ozone, stratospheric temperature,or UV radiation, a model of the followingform has frequently been found to be useful. Let Yt be our measurementon the geophysicalvariable of interest,either on its original scale or after a transformation(most often the log transformation,e.g., log UV). Then the linear trend model assumesYt = I• + St + roXt + N t, where /• is a constant term,X t = t/12 represents the lineartrendfunction,and rois the magnitudeof the trendperyear.Stis a seasonal component which
impactof suddeninterventionson the abilityto detecttrends. Suchlevel shiftshave been observedin empiricaltrend studies for a varietyof data sets[e.g.,Reinselet al., 1994;Weatherhead et al., 1997]. The basic statisticaltheory used in this paper is
canoftenberepresented asSt -- 5;•4-•[/3•,jsin(2•rjt/12)+
2.
variance• of thewhitenoiseprocess {et} in thismodelby o-} = Var (Nt)= 0-•2/(1- (•)2).
/32,jcos(2•rjt/12)]. Whiletheseasonal component is essen-
tial in practicalmodelingof geophysical time series,estimation of thiscomponentdoesnot havemuchimpacton the statistical coveredin more detail in standardtextbookson regressionand propertiesof the estimatesof the other terms in the model. time seriesanalysis[e.g.,Box et al., 1994]; however,here ap- Hence the seasonalcomponentwill not be included in our plicationsand examplesare presentedwith specialemphasis subsequentstatisticalderivations,for convenience,althoughit on the case of linear trend estimation with examplesfrom will be includedin the analysesfor the empiricaldata applicaenvironmental data sets. tions. Therefore to investigatethe effects of magnitude and The detectabilityof a trend can be summarizedin several autocorrelation of noise on trend estimates, we consider a ways.Two commonways are through the precisionof a trend simpletrend model of the form estimate, as measuredby its standarddeviation, and through Yt-- !x + roXt+ Nt, t = 1,''', T. (1) the numberof yearsof data requiredto detecta trend of given magnitude using the trend estimate. Which summaryis emFor the unexplainedportion of the data the noise Nt is ployedis determinedby what data are availableor are being assumedto be autoregressive of the orderof 1 [AR(1)]; that is, plannedfor. If the data havealreadybeencollected,generally Nt = qbNt-1 + et, where the et are independent random the former is more useful, while if an experiment is being variables withmeanzeroandcommon varianceo-•.It will also planned or a monitoring project started, the latter is often be assumedthat -1 < 4>< 1, so the noiseprocess{Nt} is more useful. In each of the following sections,the statistical stationary.This model allowsfor the noiseto be (auto)corremodelingapproachwill be introduced,and the precisionof a lated among successive measurements, such that 4> = trend estimate will be examined. Then the number of years Corr (Nt, Nt_•). The autocorrelation,which is typicallyposrequired to detect a giventrend will be discussed, followedby itive, may be the result of various natural factorswhich give exampleswith actual data sets to illustrate the precisionof rise to somewhatsmoothlyvarying changesin N t over time. trend estimates.While exampleswill be suppliedusingatmo- Suchnatural factorsmay not alwaysbe known or measurable, sphericdata, the resultsare generaland are applicableto data and the lagged value N t_ 1 in the model can sometimesbe from a variety of sources.Details of the statisticalapproaches regarded as a proxy to representthese natural factorswhich will be more fully outlined in the Appendix. dynamicallyinfluence the current noise value N t. It is also noted that the varianceof the noiseN t is directlyrelated to the Basic Trend
Evaluation:
Effects
of autocorrelation and variability on trend
estimation
and detection
Changesin environmentalvariables are often statistically modeled as being a smooth linear change. While this may never be the case, it allows a simple approximationof the directionand magnitudeof the changesin the data, which may be adequatefor many practicalpurposes.More importantly, resultsfrom the linear trend models are commonlyused and are familiar to scientistsand policymakers;thusit is important to examine various issues under such models. The linear trend
in environmentaldata is often expressedin percentchangeper decade. There is often a great deal of variation on different timescalesoccurringwithin the data in addition to the linear trend. For independent(uncorrelated)time seriesdata the length of time to detect the trend is determined by the magnitude of variation of the noise. However, environmental data
are often autocorrelated. For instance, higher than normal
2.2.
Trend
Estimation
A common situation of trend evaluation is when a given number of yearsof data have been collectedand one wantsto estimate a trend based on the existingdata. For such a situation, one wants to estimatethe trend as well as the precision
(standarddeviation)of the trend estimate.Let •o denotethe generalizedleast squares(GLS) estimatorof the trend ro in model (1), and let o-co = s.d. (•o) denote the corresponding standarddeviation of •o.The exact form of o'co is derived in the
Appendix, with•o = Var(•o)givenin (A5). In theAppendix it is alsoshownthat a usefuland quite accurateapproximationis
O'& • (1--(•))/•3/2--/•3/2 __(•) where n = T/12 denotesthe number of years of data.
(2)
WEATHERHEAD
ET AL.: TREND DETECTION
As seenfrom approximation(2), the precisionof & is directly a function of the magnitudeof variation and the autocorrelationof the noise,aswell asthe fixed numberof yearsof available data. Specifically, data with larger variance and higher positiveautocorrelationof the noise decreasethe precision(i.e., increases.d. (&)) of the trend estimate.Large o-N meanslargenoisein the data, makingany signalor trend more difficult to detect, the standard "signal-to-noise"issue. Increasesin positive autocorrelation(rk) tend to increasethe length of trend-like segmentsin the data and thus make the
17,151
no autocorrelation(assuminga commonvalue of the noise
varianceo-2•).For numerical illustration, suppose thatoneis interestedin beingableto detecta trendof magnitude5% per decade(•o0 = 0.5% per year). Table 1 (top) showsthe number of yearsof data that are neededto detectsucha trend given datawith variousvaluesof boththe magnitudeof variationo-N and autocorrelationrk.The resultspresentedshowa varietyof features.First, the number of yearsto detect a giventrend is
stronglyinfluencedby both the autocorrelationand the variance of noise in the data. For example, a month-to-month estimation of the real trend more difficult. variabilityof 10% (o-N = 10) and moderateautocorrelationof Consequently,failure in takinginto accountthe propervari- rk = 0.5 resultsin a situationrequiring 23.4 years of data to ance and autocorrelation of the noise can lead to a false level detecta 5% per decadetrend. For the samevalue of o-N but of precisionin the trend estimate,in that the assumedstandard with zero autocorrelation,4>- 0, the requirednumberof years deviationof the trend estimatewill substantiallyunderstatethe will be 16.3,a differenceof about7 yearscomparedto the case actualuncertainty.For instance,a trend result maybe obtained rk = 0.5. This result also showsthe effect of ignoring the for typicalautocorrelatedenvironmentaltime seriesdata using presenceof autocorrelationin the data (e.g., assume4>= 0 a statisticalmodel that ignoresautocorrelation,for example, when in fact rk = 0.5). There would then be a false senseof assuming rk= 0 in model(1). The resultobtainedmightfalsely confidencein the ability to detect a trend; that is, the number indicate a statisticallysignificanttrend at the 95% confidence of yearsto detecta trend would be underestimatedby 7 years. level, whereas if the appropriate statisticalmodel were used Note that the resultsin Table 1 (top) are alsopreseptedin whichdoesnot ignorethe autocorrelation,the actualprecision termsof the ratio O-•v/C00 and can be usedfor generalvaluesof of the trend estimatemight be found to be substantiallyless. trend (c00).For illustration,if it is desiredto detecta trend of For any givenvalue of 4>,we seefrom the approximationin magnitude •o0 = 0.25% per year, and the underlyingnoise (2) that o-&/o-•decreasesas n increasesat a rate essentially processhas O-•v= 3.0%, then we can obtain the required proportional to/./-3/2. The standard deviationof & increases numberof yearsfrom the row in Table 1 (top) with O-•v/C00 = as rk increases,and with a fixed value of o-•, the standard 12 (e.g., n* = 16.6 when rk = 0.5). deviation of & for the moderate autocorrelationof rk = 0.5 is For further numericalcomparison,Table 1 (bottom) shows approximatelydoubledrelative to its value when rk = 0. How- the number of years n* of data required for detection of a ever, it mustbe noted that AR(1) noiseprocesses{N,} with trend of magnitude1% per decade;that is, c00= 0.1% per year. the sameo-• - s.d. (et) and different rkvalueshave different Even for small autocorrelationand low variability, it would variances,since O'N 2 = Var(Nt) = o-•2/(1- qb2).Hencein often require a prohibitivelylong period of time to detectsuch termsof o-N the standarddeviationof & canbe written as o-co = a smalltrend. The longtimesrequiredfor trend detectionarise O'NV' (1 + rk)/(1 -- rk)n-3/2. Thusif o-N, ratherthano-•,is becauseof the random variability in the data sets.This varifixed, then the standarddeviationof & when 4>= 0.5 is in fact ability is itself of great scientificinterest. For many environapproximately V'(1 + 0.5)/(1- 0.5) = X/-J• 1.73timesthat mental parameters,for example,heightof sealevel, precipitawhen 4>= 0. tion, or UV radiation, a year of unusuallyhigh levelsonce in a while mightbe of greaterrelevancethan a smallupwardtrend 2.3.
Trend
over many years.
Detection
The issueof key interest in trend detection is to determine the number of years of data required to detect a trend of a certain magnitude.We shall adopt the commonlyused decision rule that a real trend is indicated,at the 5% significance
2.4.
Uncertainties
Due To Unknown
Variance
and Autocorrelation
The approximationfor n * givenin (3) not only assumes the linear model (1) is adequatebut also that we know the true level or 95% confidence level,when Ico/l > 2. It is then valuesof O-•vand 4>-In practice, even if we assumethe adeestablishedby Tiao et al. [1990]that the numberof yearsn * of quacyof model (1), we still do not know O-•vor rk,and thuswe data required to detect a real trend of specifiedmagnitude haveto estimatethem from availabledata (and prior informa-
l01- I001,with probability 0.90,is 3.3o-•
I00( -
]2/3 :
tion);thatis,wewillusetheestimated values (6-•and•b)in
,
(3)
placeof o-•and 4>,respectively, whenwe apply(3), resultingin an estimatedn*, •* Given that n* stronglydependson o-•
and rk,it is naturalto worryaboutthe uncertain•ty in •i*,
wherewe assume101< 1. We note here that this formula particularlyif there is large uncertaintyin 6-• or in 4>.In pracrequires O-•vand •o0 be on the same scale;that is, if •o0 is in percentage,then O-•vmust be convertedto percentageas well (e.g., dividingthe original O-•vby a mean). Also, the value 3.3 in (3) wouldbe larger(smaller)if we requirea higher(lower) degree of certainty than 90%. From this result we can see the importance of both the magnitudeof variation and the autocorrelationof the noiseon the detectionof trend. In general,from (3) we find that the presenceof positive autocorrelation(rk > 0) increasesthe required number of years for trend detectionby the factor of
[(1 + 4>)/(1- qb)] •/3,approximately, compared to thecaseof
tice, only a few yearsof data are necessaryto estimateo-•fairly well, and thusit is often acceptableto ignorethe uncertaintyin &• (i.e., we can treat o-• = 6-•); significantlylonger time is
needed toadequately estimate 4>.Consequently, it isimportan•t to assessthe uncertaintyin •i* due to the uncertaintyin rk, consideringthat the amount of data availablefor estimatingrk at the planning stage (i.e., when we need to estimate the numberof years) is typicallynot large.
Thefollowing method pr•ovides anapproximate 95%confidenceintervalfor n * when rkis the leastsquaresestimateof qb based on M months of data. The derivation is given in the
17,152
WEATHERHEAD
ET AL.: TREND DETECTION
Table 1. Number of Years of Monthly Data Needed to Detect, With Probability0.90, a Trend
of 5% Per Decade
and 1% Per Decade
at a 95% Confidence
Level for Selected
Values of Autocorrelation(qb) and StandardDeviation (rrN) of the Noise Value of 4>
O'•v
(O'er/W0)
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.25 0.5 1 2 4 6 8 10 12 15 20
Number of Yearsto Detecta Trendof 5% Per Decade (•oo = 0.5% per year) (0.5) 1.4 1.5 1.6 1.7 1.8 2.0 2.1 2.3 (1) 2.2 2.4 2.5 2.7 2.8 3.0 3.3 3.6 (2) 3.5 3.7 4.0 4.3 4.6 4.9 5.3 5.9 (4) 5.6 6.0 6.4 6.8 7.3 7.9 8.6 9.6 (8) 8.9 9.5 10.1 10.8 11.6 12.6 13.8 15.4 (12) 11.6 12.4 13.3 14.2 15.3 16.6 18.2 20.3 (16) 14.1 15.0 16.1 17.2 18.6 20.1 22.1 24.7 (20) 16.3 17.4 18.7 20.0 21.6 23.4 25.7 28.7 (24) 18.4 19.7 21.1 22.6 24.4 26.4 29.0 32.5 (30) 21.4 22.9 24.5 26.2 28.3 30.7 33.7 37.8 (40) 25.9 27.7 29.6 31.8 34.3 37.2 40.9 45.8
4.0 6.7 11.0 17.8 23.5 28.6 33.3 37.7 43.8 53.3
0.25 0.5 1 2 4 6 8 10 12 15 20
Number of Yearsto Detecta Trendof 1% Per Decade (•oo = 0.1% per year) (2.5) 4.1 4.3 4.6 4.9 5.3 5.7 6.2 6.9 (5) 6.5 6.9 7.4 7.9 8.5 9.2 10.0 11.2 (10) 10.3 11.0 11.7 12.6 13.5 14.7 16.1 18.0 (20) 16.3 17.4 18.7 20.0 21.6 23.4 25.7 28.7 (40) 25.9 27.7 29.6 31.8 34.3 37.2 40.9 45.8 (60) 34.0 36.3 38.8 41.7 44.9 48.8 53.7 60.2 (80) 41.2 44.0 47.1 50.5 54.5 59.2 65.1 73.0 (100) 47.8 51.0 54.6 58.6 63.2 68.7 75.6 84.7 (120) 53.9 57.6 61.7 66.2 71.4 77.6 85.4 95.8 (150) 62.6 66.9 71.6 76.8 82.9 90.1 99.1 111.2 (200) 75.8 81.0 86.7 93.1 100.4 109.2 120.1 134.8
7.8 12.8 20.7 33.3 53.3 70.0 84.9 98.7 111.5 129.5 157.0
The standarddeviation o-N of the noiseis expressedin percentvariability in the month-to-monthdata. Table valuescan also be usedfor trend detectionfor generalvalue of •o0,in terms of the ratio crN/•O o.
Appendix(sectionA2). Firstwe calculatean estimateduncer- 2.5. Applications tainty factor
The precisionof a long-termlinear trend estimatehasbeen shownin sections2.1 and 2.2 to be highly dependenton the
B=3x/• _•'
variance(• or rr2N)and autocorrelation (4>)of the noise.
(4) Realisticestimatesof theseparameterscan allow one to esti-
mate the number of yearsnecessaryto detecta trend of given given by(• *e-B, h*eB),where g* iscalculated according to magnitude,or the magnitudeof trend that can be detectedby (3) using4>and &• for qband cry,respectively.Note that when a fixednumberof yearsof data. Knowledgeof theseparameter B _< 0.2, the intervalis essentiallygivenby h*(1 _+B), thus values(rrN and 4>)overa varietyof differentlocationscanhelp B is the percentageof uncertaintyrelativeto the point estimate one determine the most likely locationsfrom which a trend g *. For largerB theinterval(• * e- B, •. eB) is notsymmetric could first be detected.Suchan analysiswill allow reasonable estimatesof the time necessaryto detect a giventrend. about • * but rather skewedtoward large values. To illustrate, UV data from 14 Robertson-Berger(RB) To illustratethe calculationof this procedure,suppose6-• = Then the 95% confidenceinterval for the number of years is
3.1%,co o = 0.3%peryearand•b= 0.32.Thesevaluesare
meter
stations within
the United
States were examined
to as-
actually from a total ozone data set collected at Tateno, as sessthe autocorrelation4>and the magnitudeof variabilityrrN reported by Tiao et al. [1990], who gave •* • 14 years;the for trend detectionof UV radiation after seasonalvariability et al., 1997].Table 2 valuefrom(3) is 13.6.Nowsuppose thevalue•b-- 0.32was hasbeen accountedfor [seeWeatherhead liststhe 14 stationsand givesthe estimatesof 4>and rrN and the computed on the basis of 2 years of data; that is, M = 24. Then the B givenby (4) is 0.38, and thusan approximate95% correspondingnumber of years required to detect a trend of 5% per decade.The resultsare also graphicallyillustratedin intervalfor thenumberof yearsis (13.6e-ø'38,13.6eø'38)• (9 years,20 years),whichis a rather wide interval.If we have Figure 1, which showsthe estimatesof 4>and rrN alongwith 5 yearsof data to estimate4>(i.e., M = 60), then B = 0.24, contour curvesindicatingthe number of yearsof data needed and the intervalbecomes(10.8, 17.3), or approximately14 _+3 to detecta 5% per decadetrend. The plot showsthat the range years.Having 5 yearsof data to estimate4>duringthe planning of autocorrelationand variability observedat these stations stageis likely to be unusualin practice. AssumingM varies translatesinto a rather wide range of the number of years to from 24 to 60 months and 4>varies from 0 to 0.5, the relative detect a giventrend. It shouldbe noted that the valuesof rrN uncertainty from the interval estimate varies approximately and 4>havea latitudinaldependencewith lower rrN and higher from 20 to 60%, and thus it is clear that the uncertaintyin 4>observednearer the equator(seeTable 2). estimating4>cannotbe ignoredwhen estimatingthe numberof While suchan analysiscanbe helpful in identifyinglocations years. likely to detect a given trend earliest,it shouldbe noted that
WEATHERHEAD
ET AL.: TREND DETECTION
17,153
Table2. Estimated Standard Deviation (&N)andAutocorrelation (•b)of Noisein UV Data From 14 LocationsAlongWith the EstimatedNumberof Years (•*) to DetectTrends Reference
&N
Latitude
Site
Trend,
•b
(in %)
Seattle, WA Bismarck, ND
47.5øN 46.8øN
0.12 0.24
13.0 14.9
% perdecade
g * for
g * for 5%
Reference
Trend
3.1 2.6
29 37
Trend,
perdecade 21 24
Minneapolis,MN
44.9øN
0.13
11.0
3.0
27
19
Detroit, MI Des Moines, IA Salt Lake City, UT
42.4øN 41.6øN 40.8øN
0.12 0.29 -0.15
11.6 13.5 10.5
3.2 2.9 2.3
26 35 25
20 24 15
Philadelphia,PA
39.9øN
0.05
12.3
2.9
28
19
Oakland, CA
37.7øN
0.21
9.3
2.3
30
18
Albuquerque,NM
35.1øN
0.31
7.7
2.1
30
17
Fort Worth, TX Tucson, AZ E1 Paso, TX Tallahassee, FL Mauna Loa, HI
32.8øN 32.3øN 31.8øN 30.4øN 19.5øN
0.23 0.27 0.18 0.04 0.36
12.4 7.0 7.8 9.8 12.0
1.7 1.9 1.7 1.8 0.6
45 30 33 33 97
22 16 16 17 24
The uncertaintiesfor the estimatednumber of years presentedin this table are all of the order of _+10-15% becausethe number of yearsof data here range between 10 and 18 years.
for anygivenvariable,in thiscaseUV radiation,the magnitude of the expectedtrend mayvaryby locationaswell. Thus a more pointed analysiscomparingsize of expectedtrendswith the ability of a site to detect that trend is of further value. Often a realistic expectedtrend may not exist, but a reference trend may be establishedusing the best available information. In practice, the estimatednumber of years to detect a trend is highly sensitiveto the magnitudeof the reference trend. In practice,a sensitivitystudycan show the exact dependence. Columns6 and 7 of Table 2 illustratethe importanceof estab-
lishing expectationsfor trend detectionin evaluatingvarious sites for their ability to detect trends. On the basis of past ozonetrends,UV radiation is not expectedto changeas much in the lower-latitude sites as in the higher-latitude sites. In particular,usingpast ozone data from the total ozonemapping spectrometeron the Nimbus 7 satellite as a reference and assumingthat UV radiation will increaseapproximately1.9% for each 1% decrease in ozone, we can establish reference
trendsfor UV radiationas givenin Table 2 (column5). Calculatingthe numbersof yearsto detect thesereferencetrends, shown in column 6 of Table 2, reveals that data from Detroit
Measured
Ultraviolet
and Minneapolis are likely to identify their reference trends severalyearsearlier than Tallahassee,E1 Paso,or Tucson.This is in contrastto the resultspresentedin column 7 of Table 2
Radiation
which show that data from sites like Tallahassee, E1 Paso, and
Tucsonwould be likely to detect a commontrend of 5% per decade several years sooner than from sites like Detroit or Minneapolis.A more completecomparisonwould includethe
(•
estimate
(•
0
10
15
20
yrs
yrs
yrs
10
yrs
20
yrs
30
40
of errors in these numbers.
In addition to the comparisonof detectabilityat different locations,there is another comparisonthat is of great interest whichthis type of analysiscan offer. Anthropogenicinfluences are predictedto create changesin a variety of variables.The variablewith the largestexpectedchangemay be identifiedas a key variable to monitor with respectto globalchange.However, different environmentalvariableshave characteristically different noise associatedwith them, implyingthat for particular variables,despite their low expectedsignal, it may be possibleto detecta long-termchangerelativelyearly.Figure 2 showsestimatesof o-N in percentof meanvalue and qbfor data on sixatmosphericvariablesafter seasonalvariabilityhasbeen accountedfor, with the data for the first four plots obtained from the same 14 locationsin the United States,as given in Table
2. The surface ozone
data are from
13 locations
inter-
nationally,as listed in Table 3. The total ozone data are from satellite measurementsover the same locations,as presented Figure 1. Magnitude(O'N) and autocorrelation(qb)from 14 in Figure 1 and Table 2. Informationon the sourcesof the data UV monitoring sites within the United States. Lines show combinations of o-N and qbrequiredto detecta trend of 5% per maybe foundin the Appendix(sectionA1). The linesin Figure decade in the stated number of years. This plot showsthat 2 indicatethe numberof yearsof data neededto detecta trend somelocationsappearto be better suitedthan other locations of 5% per decade.Figure 2 showsthat some environmental for detectinga given trend. parameters appear to be inherently more variable or more N
17,154
WEATHERHEAD ET AL.: TREND DETECTION
Trend Detectabilityfor a Variety of Parameters Solar Radiation
o
Solar Radiation
- Horizontal
- Direct
_
10
15
20
30
40
10
yrs
yrs
yrs
yrs yrs
yrs yrs i
i
i
i
0
10
20
30
40
15
4O
20
yrs
yrs
i
i
i
i
i
0
10
20
30
40
N
N
Surface Pressure
Relative Humidity
c
10
15
yrs yrs i
i
i
0
10
20
20
30
40
yrs
yrs
yrs
i
l
i
i
i
30
40
0
10
20
Surface
30
40
30
40
N
N
Total Ozone
Ozone
ß ß
(5
10 I
;
yrs
yrs
yrs
yrs
i
0
yrs
yrs
yrs
y
,
10
20 N
30
40
0
10
20 N
Figure 2. Magnitude(O'N)andautocorrelation (45)from a varietyof differentmonitoringstations. The first four plotsshowdatafrom the same14 locationswith the followingkey for locations:A, Caribou,Maine; B, Bismarck,North Dakota; C, Great Falls, Montana; D, Eugene,Oregon;E, Lander, Wyoming;F, Madison, Wisconsin;G, Omaha,Nebraska;H, Sterling,Virginia; I, Albuquerque,New Mexico;J, Phoenix,Arizona; K, E1 Paso,Texas;L, Montgomery,Alabama;M, Tallahassee,Florida;N, Brownsville,Texas.(Letter order representsdecreasing latitude.)The plot of surfaceozonerepresents datafrom 13 locationsinternationally. The plot of total ozonerepresents datafrom satellitewhichcorrespond to the 14locationsdisplayed in Figure 1. In all sixplots,linesshowcombinations of trN and & requiredto detecta trendof 5% per decadein the statednumberof years.These plots showthat somevariablesappearto be inherentlymore difficultthan others for detectinga fixed trend.
autocorrelated than others. The variables examined for this tions showsthat while a site, suchas Tallahassee,Florida, may study indicate that some are capable of showingsignificant be a preferred site for detectinga trend in one parameter changessoonerthan others.Examinationof the firstfour plots (relativehumidity),it maybe a lessdesirablesitefor detecting which showa number of parametersfrom the same 14 loca- a trend in anotherparameter(solar radiation).This analysis
WEATHERHEAD
ET AL.: TREND DETECTION
17,155
Table3. Estimated Standard Deviation (d-N)andAutocorrelation (/b)of Noisein Surface OzoneData from 13 LocationsAlongWith the EstimatedNumberof Years (•*) to Detect Trendsand an Approximate95% Interval for the Number of Years (in Parenthesis) h*
to Detect
5% per Decade Trend
Site Barrow, Alaska
M
•b
%
270
0.56
39 78
0.64 0.41
204 53 83
Iceland Ireland
Zugspitze Niwot Bermuda
h*
to Detect
1
ppb per Decade
and Its
and Its
ppb
Uncertainty
Uncertainty
22.5
5.77
11.5 10.7
4.35 3.87
0.61
10.6
4.17
0.74 0.65
9.2 20.7
4.19 7.06
-42 (36, 50) 29 (18, 47) 23 (17, 29) 27 (22, 33) 28 (17, 46) 44 (32, 61) 32 (23, 43) 42 (34, 51) 27 (19, 36) 56 (45, 69) 27 (18, 39) 12(9,15) 45 (35, 56)
50 (43, 59) 45 (28, 71) 34 (26, 43) 43 (35, 52) 49 (30, 80) 63 (45, 87) 55 (40, 75) 63 (52, 77) 27 (20, 36) 44 (36, 54) 27 (19, 39) 14(11,17) 58 (46, 73)
Canary
104
0.68
12.1
5.46
Mauna Loa Barbados Samoa
259 53 214
0.67 0.43 0.65
18.5 13.3 29.4
6.88 3.68 4.14
Cape Point Cape Grim
121 154
0.79 0.56
South Pole
243
0.75
7.6 3.4 17.3
1.54 0.85 5.16
Also listedis the length of the data (in month,M) usedfor estimation.
showsthat in termsof detectinga statisticallysignificanttrend, there is no singlelocationthat is bestfor all parameters.Again, the expectedchangefor a variableneedsto be comparedto the likelihood of detectinga given changein a variable to determine the variablesmostlikely to be usefulfor the detectionof a trend. In the absenceof expectedor reference trend values, in Figure 2 all sitesare beingcomparedon the basisof detecting a 5% per decadetrend. The presentationin this sectionassumesthat one is looking for a trend of a givenpercentper decade.Often, the expected changeis not a percentagechangebut may be representedas an absolutechange,suchas an increaseof surfaceozone of 1 ppb per decade.Table 3 showssome analysisof the surface ozone data also presentedin Figure 2, where the data were shownin terms of detectinga 5% per decadetrend. Scientifically,it may be more appropriateto look for a fixed changein concentration
of surface
ozone.
Table
3 therefore
shows the
numberof yearsto detect a 5% per decadetrend aswell asthe number of years to detect a 1 ppb change in surface ozone. Note that the data from Ireland appearto be more appropriate for detectinga 5% per decadetrend than the data from Cape Point; however,the data from Cape Point are more appropriate for detectinga trend of 1 ppb per decade than the data from Ireland. Realistically,the uncertaintyin the estimatefor the number of years for these two sites are not statistically distinguishable.However, this example is only illustrative in nature.
These analysesshowthe importanceof estimatinghow long it will take to detect trends at different
locations
and in differ-
ent environmentalparameters.Proper analysisof trend detectabilityfor different locationscan help focusexistingand future monitoring activities.Choosingsites appropriatelycan allow scientificquestionsrelevant to the health of the environment to be answered sooner. Efficient monitoring can also save considerablyin termsof coststo monitor and delaysassociated with waiting for reliable trend results.
for decades.However, often unforeseenproblemsarisewhich disruptthe measurementsin someway that cannotbe quantified with any certainty. Instruments may break or may be modified;sitesmay changelocation, maintenanceand calibration proceduresmay change,or external eventssuchas volcanic eruptionsmay occur.A usefulstatisticalmodel for manyof theseeffectsis obtainedby assumingthat at a specifiedpoint in time the data experiencea permanent level shift of unknown magnitude [Box and Tiao, 1975]. If the data are analyzed in logarithmicform, sucha discretelevel shift in the logarithmsof the original data correspondsto a scalefactor change in the original data that could represent,for instance,the changein sensitivityof an instrument.To ignore sucheventscan result in artificial trends in the data that are not representativeof the environmentalvariable being studied. If the magnitude and time of the shift are known, then the data can be adjusted before being analyzedfor trends. However, often the time of the level shift is knownbut not the magnitude.The presenceof such a level shift will result in an increase
in the variance
of a
trend estimatethat properly accountsfor the estimatedshift and hencelengthensthe time necessaryto detect a giventrend. The impact of the level shift on trend detectionstronglydependson its relative location in time in the data set. 3.1.
Basic Statistical Modeling
To investigate,statistically,the effectsof interventionlevel shift, in addition to autocorrelation, on trend estimates, we considera trend model of the following form:
Y, = /• + coX,+ 8U, + N,,
t = 1, ß-., T,
(5)
where foxt representsa linear trend beginningat time t = 0, as before, and BUt is a mean level shift term used to account for the possibleinterventionto the data at the specifiedtime t = T O (0 < T O < T); that is,
->To' 0, tt< To
Ut= 1, 3.
Trend
Evaluation
With
Interventions
The noiseseriesN t is modeledas an autoregressive [AR(1)]
Data collection for the detection of trends often requires process,Nt = &Nt_ •. + et, as in model (1). We are particumaintenanceof stableinstrumentsin representativelocations larly interestedin the uncertaintyor varianceof the estimateof
17,156
WEATHERHEAD
ET AL.: TREND DETECTION
If one ignoresa potentiallevel shiftin the statisticalmodel,the resultingtrend estimatewill be biased.The exactform of the bias is presentedin the Appendix.As an example,in the case of no autocorrelation with 4>= 0, it reduces to a bias in the extent on whether the time of the intervention is known and derivedtrend proportionalto the magnitudeof the level shift: bias in do • 6,(1 - ,)/n 15.In practice,level shifts of whether any additionalinformationis availableregardingthe magnitudeof the level shift. For the presentwe considerthe moderate magnitude can overwhelm and obscure the true et al. [1997]and casewhere the time of interventiont - T O in model (5) is trend beingsought,as shownby Weatherhead known. In somecases,althoughthe existenceof a (possible) Krzy•cin [1996] in the analysisof UV data. In addition, in interventionwithin a certaintime interval maybe expected,the practicethe estimateof the white noisevariance• will be exact time T O of the interventionis not known.A variety of biasedupwardwhen the level shift term is not includedin the model. statisticalwork has been done to investigatethis problem, includingestimationof the interventiontime T O and testing the null hypothesisof no intervention.Generally,it is found 3.3. Extension When There is Additional Information that propertiesof estimatorsof the trend and shiftcoefficients In somesituations,additional(prior) informationmay be (to and/5) in model(5) are not affectedmuchwhen T Ohasto availableon the magnitudeof the mean level shift/5 in model be estimated. The case of additional information about the (5). In practice,when the level shift couldbe causedby instrument-relatedproblems,for example,replacementor recalibramagnitudeof the level shiftwill be discussed in section3.3. It is importantto emphasizethat model(5) is not universally tion of instrument,this additionalinformationmay comefrom applicablefor all typesof interventions.Long-termdrifts in some comparativemeasurementsof the instrument againsta instrumentationwhich are periodicallyadjustedand volcanic standard instrument, or other information available from calinterventionswhere the effect of the eruptiondecayswith time ibration procedures,a measurementcampaignusing two inare exampleswheremodel(5) wouldnot be appropriate.How- strumentsfor some overlappingtime period surroundingthe ever, local changesin sites,maintenanceprocedures,calibra- time of intervention. To formalize the nature of information, tion techniques,and instrumentationmay be adequatelyde- statistically,we assumethat for the magnitudeof the level shift scribed by the model considered here. The use of an (/5) we have a (prior) estimate/50with specifiedvariance(a
the trend in model (5), and the effect of occurrenceof the interventionmean level shift term (as well as the autocorrelated AR error) on the varianceof the estimatedtrend. The statisticaltreatment of level shiftsdependsto a certain
inappropriatemodel is likely to result in biasedestimatesof trends as well as other parameters. 3.2.
Trend
measure of the precision of priorestimate)•. Typically, the prior estimate/50and its variance• will be obtainedfrom independentanalysisof the comparativedata, for example, simultaneous measurements from two different instruments,
Estimation
mentioned
above. We also assume that no other information
is
To describethe impactof interventionson trend estimation, let dodenotethe GLS estimatorof the trend toin (5), and let available or incorporated into the analysisconcerning the trco= s.d. (do)denotethe standarddeviationof do.The exact other regressionparameters/xand toin model (5). We denotethe ratio• = •/• asa relativemeasure of the form of trcois derivedin the Appendix,and the expressionfor extentto whichwe know the magnitudeof the level shift.Then, thevariance (•0,)of doisgivenby(A4).A usefulapproximation as shownfrom (A6) in the Appendix, the variance of the for the standarddeviation of the trend estimateis given by correspondingtrend estimatedois given by the expressionin (A4) but with h 6 replacedby h 6 q- /(' From this we find that the variance
O'dø•(1 -- 4))Ft 3/2 [1 -- 3T(1-- T)]1/2
of the trend estimate
can be decreased
substan-
tially asthe amountof additionalinformationincreases(i.e., • increases).To illustrate,considera seriesof length n - 10 = n -• [1- 3,(1- ,)],/2, (6) years with the intervention halfway through the time series (•- = 0.5), and moderateautocorrelationof 4>= 0.5. When where •- - (T o - 1)/T is the fraction of data before the • = 0, 1, or 9, the standarddeviation of the trend estimate trco intervention.This approximationis accuratefor smallervalues is about 1.8, 1.6, or 1.2, respectively,times that when no interof (kbut tendsto overstatethe true value of trco as (kgetslarger; ventionis present(•- = 0). Thuswhenadditionalinformationis also,it becomesmore accurateas the number of yearsof data available correspondingto a value of • - 9, the standard n getslarger. When there is no intervention,•- = 0, the result deviation of the trend estimate is substantiallyreduced over in (6) reducesto that givenin (2). The effect of a level shift the caseof no additionalinformation(• = 0) and yieldsvariintervention on the standard deviation of the trend estimate ability not too much greaterthan the casewhere no interven(trco) is represented bythefactor1/[1- 3,(1 - ,)]•/2,regard- tion is involved,equivalently,where the magnitudeof the inlessof the value of (k.From (6) we seethat the uncertaintyof terventionis knownwith certainty. a trend estimateis greatestwhen the interventionoccurshalfway through the collectionof data (•- -- 1/2); in this case, 3.4. Trend Detection
1/[1- 3,(1 - ,)]•/2 = 2, sothestandard deviation trco is2 times the standard
deviation
obtained
in the no-intervention
case
from (2). When •- = 0.25 or 0.75, the standarddeviation&,, is about 1.5 times the value when no interventionis present. While the resultsaboveshowthat includinga level shiftterm in the statisticalmodel usedwill increasethe uncertaintyof the trend estimate, more seriousconsequencesoccur if one does not includea level shift term when a level shifthastaken place.
Similar to the no-intervention case considered in section 2.3,
the number of yearsn* of data required to detect a trend of magnitudeto0 can be determinedunder the level shift intervention model (5). This model assumesthat no additional information on the magnitudeof the level shift is available. Using the approximationfor the standard deviation of the trend estimategivenin (5), a rough but convenientapproximation for the number of years for detection is obtained as
WEATHERHEAD
3.3•r•
ET AL.: TREND
I2/3
n*= w0(1- 4))[1 - 3,(1- ,)]•/2
DETECTION
17,157
estimation
results are shown in Table
4. We see that in each
case there
is a substantial
in trend
tween the no-level
difference
estimates
be-
shift and the level shift models.
number of years for trend detection roughly by a factor of
As discussed in sections3.1 and 3.2, generally,the inclusion of a level shift term has the effect of increasingthe standard deviationof the trend estimate,when other thingsare equal. However, in practice, if a substantiallevel shift intervention existsbut is not accountedfor in the statisticalmodeling,then the error varianceis undulyinflated.Thusinclusionof the level shiftterm in suchprominentcasesof interventionwill typically
1/[1- 3,(1 - Q]l/3.For example, for theworstcase,i.e.,, =
result in a decrease
1/2, this approximationgivesa hctor of 1.59. Using the exactexpressionfor the standarddeviationof the trend estimategiven in (A4), the number of years (n*) requiredto detecta 5% per decadetrend (wo = 0.5% per year), for the case of an inte•ention hal•ay through the data set (, = 1/2) are obtainedfor different values of a• and &. By comparing these values with the correspondingvalues from Table i (top), we find that the ratio of the numberof yearsfor
timate6-•.It mayalsoleadto a reduced estimate •bof the
Comparisonof thisconvenientapproximationwith the approximation (3) for the no-inte•ention casesuggests that the presence of a level shift inte•ention
trend
detection
no-inte•ention
in the model
increases
to that under
the
under
the level shift model
model
is about 1.55 to 1.59 for smaller values
the
in the white
noise standard
deviation
es-
AR(i) coefficient,sincesomeof the low-frequencybehaviorof the noise that was actually causedby the level shift is then explicitlymodeled in the level shift term and hence lesslow frequency autocorrelationis present in the remaining noise
term.As a resultof thesepossible reductions in 6-•and•b,in practicethe standarddeviation of a trend estimatewith a level shift term includedmay not be muchlarger (or couldevenbe smaller) than that obtainedfrom the model without the level shift term. For the results shownin Table 4, these features tend to be present to some extent. However, ideally, if the level shift feature had not been presentin the data (or couldhave been "eliminated"in some
of & (& • 0.5), andthe ratio becomesslightlysmallerfor larger valuesof &. This ratio is not ve• sensitiveto the value of a• unless& is large. Thus the presenceof an inte•ention level shifthal•ay throughthe data collectioncausesroughlya 50% justifiable way),andthedatapossessed thesame6-•and•bas increasein the numberof yearsreq•red for detectionof trends. obtained from the above level shift model estimation results, We note that the proceduregivenin section2.4 for assessing then a trend estimate from the no-level shift intervention the uncertainWin the estimatednumber of years,• *, due to model would have a substantiallysmaller standarddeviation, the uncertaintyin estimating& is applicablewhen • * is calcu- as indicated by the results of section 3.2. These "idealized" lated according to (7). This is because the extra factor standard deviation values of the trend estimate under a no-
level shift model are shownas •r•, in column8 of Table 4. 1/[1- 3,(1- ,)]1/3 does no[depend on&andthelarge-sample formulafor the varianceof & is the sameundermodels(1) and (S).
Correspondingvaluesof the estimatednumber of years•i * for detectionof trend of magnitudewo are alsogivenin Table 4 for the "idealized"
3.5.
Applications
Consider
first the characteristics
of UV
radiation
data col-
lected from 14 RB meter stationsover the period 1974-1991. In the trend analysisof theseUV data by Weatherhead et al. [1997], usingmodelssimilarto thosediscussedearlier, it was found that the estimatesof over 13 of the 14 stationsrange from about zero to 0.35 with a "Wpical" value of about 0.2. •so, in units of 100 times the logarithmic data, a "Wpical" value for the estimate of a• is about 10. The trend estimate using such data can be interpreted roughly as a percentage trend per year. Thus, for example,if we speci• that a true trend of magnitudewo = 1% per year (as may be expectedin the •ctic) is to be detected,thenwith & = 0.2 and • = 10,we find from (7) that the required number of years of data for detectionof trend of the magnitudewo = 1% per year is about n* = 11.9 for the case of no inte•ention, = 0, n* = 15.7 for the casesof,: 0.25 or 0.75, and n* = 18.9 for the case of,
= 0.5. Hence we see from this illustration
that the occur-
rence of level shift in data resultsin several additional years required for the detectionof trend. We further
illustrate
the effect of inte•ention
level shifts on
trend estimatesand their precisionusingtotal ozone data from Huancayo, 100 mbar temperature data from Hao, and UV radiationdata from •buquerque. The deseasonalized data for each time seriesare shownin Figure 3. In each case, a level shiftinte•ention seemshighlyprobable,althoughthe causeor source is not known for all three cases.For each station, trend estimates
were
obtained
no-intervention
and intervention
cases. A more
completecomparisonwould be basedon interval estimatesfor
for both a no-level
shift and a level
shift model, which also include a seasonalcomponent, and
n *, which we omit here as we have illustrated suchcalculations before.
4.
Conclusion
In this paper, the primary statisticalconsiderationsin detectinglong-term,linear trendswere presented.The two major statisticalfactorsgoverningtrend estimationand detectionare the autocorrelationand variance of the noise. Results presented show that the number of years of data required to detect a given trend is highly dependent on both of these parameters.Examination of environmentaldata showsthat these two parameterscan vary substantiallyfrom site to site [e.g.,Zerefoset al., 1997]with the resultbeingthat one sitemay require many more yearsof data to detect a given trend than another site. Additionally, some variablesshow higher autocorrelationandvariabilitythan others,implyingthat trendsare more difficult to detect in some environmental parameters than in others.
Examination of a variety of environmentaldata setsshows that level shiftsare commonin long-term monitoring projects as instrumentsor site locationsmay change,aswell as calibration or maintenancetechniques.The work in this paper shows that the occurrenceof level shiftsin the data can add significantly to the uncertainty in trend estimatesand thereby increasethe number of years necessaryto detect a trend by as much as 50%. However, the effect of level shifts can be minimized in some situations when additional
information
is avail-
17,158
WEATHERHEAD
ET AL.' TREND
DETECTION
DeseasonalizedGeophysicalData With Level Shift Interventions Deseasonalized DobsonTotal Ozone at Huancayo, Jan. 1978 - Apr. 1992
I
I
i
i
i
i
i
i
i
i
74
76
78
80
82
84
86
88
90
92
Year
DeseasonalizedRawinsonde100mb Temperature at Hao, Sep. 1971 - Dec. 1991
!
i
i
I
i
i
!
i
i
i
i
72
74
76
78
80
82
84
86
88
90
92
Year
DeseasonalizedUV (100*log) Data at Albuquerque,Jan. 1978 - Sep. 1990
i
i
i
i
i
i
i
i
i
i
72
74
76
78
80
82
84
86
88
90
92
Year
Figure 3. Deseasonalizedgeophysicaldata with level shift interventions.The plots showdata with identifiable level shifts.Level shiftsmay indicatean overallchangein the parametermeasured,a local changenot indicativeof the parameter over a larger area or a problem in instrumentation.
able in order to estimatethe magnitudeof the level shift; and the effect can essentiallybe removed if the magnitudeof the level shift is known to a high degree of certainty. Such additional information could come from an overlap or crosscalibration in instrumentation,when changesin instrumentsare made. Estimates are made to show how varying amountsof information mitigate the impact of level shiftson trend estimation. This analysiscan be applied to a variety of practical decisionssuchas to determinethe optimal numberof months for satelliteoverlap or calibrationroutines. The resultspresentedin this paper are relevantto realistic planning of monitoring sites establishedfor the detection of trends. While there are a variety of other considerations,includingthe relevanceof a site to the scientificquestionbeing asked and the logisticmaintenanceof a monitoringsite, the statisticalconsiderationspresentedin this paper may be useful for determiningreasonableexpectationsfor trend detectability
at different sites.The applicationspresentedalsoindicatethe great need of additionalexplanatorydata and focusedstudies to detect environmentalchange.Small changesof a few percent per decadewill take prohibitivelylong time periods to detect trends, often more than several decades. Further assis-
tance in the detectionof trendsmay be obtainedfrom the use of networks,rather than singlemonitoringsites.Savings,in terms of number of years to detect a trend, depend on the location and cross correlation
of the data from the different
sites.To what extentnetworksand additionalexplanatorydata can help to reducethe numberof yearsfor detectionis a topic currentlyunder investigation. Determining which environmental variables and locations are likely to allow for earlier detectablechangeswill allow currentandplannedmonitoringprogramsto be more efficient in answeringscientificquestions.Efficient detectionof trends through focusedmonitoring activitiescan allow for scientific
WEATHERHEAD Table
4.
Illustration
of Effect
ET AL.' TREND
of Intervention
Level
DETECTION
Shift on Trend
17,159 Estimate
and Its
Precision,Basedon GeophysicalData From Three Stations
Station Huancayo(ozone in DU) Model Model
Latitude
0.715 0.527
Model Model
to
b'co
3.816 3.574
- 11.260
- 1.068 -0.040
0.237 0.228
0.137
•* = 12
g* = 8
18.1øS 0.637 0.552
with shift
(Shift at 4/76; r = 0.23) Albuquerque(UV in %)
•
(1.976)
Model without shift Model
b-•
12.1øS
without shift with shift
(Shift at 11/82;r = 0.34) Hao (temperaturein øC)
•b
1.563 1.513
-3.132
-0.164 -0.001
(0.705)
0.046 0.052
0.036
n* = 29
n* = 23
35.1øN
without shift with shift
0.277 0.193
6.993 6.775
-9.328
(2.546)
(Shift at 9/86; r = 0.69)
-0.854 0.098
0.212 0.320
0.189
g* = 18
•* = 13
Estimatesandstandarddeviationsin Dobsonunits(DU) for total ozone,øCfor temperature,and % for
UV radiation. Notation: •b,estimated autocorrelation; &•estimated standard deviation ofthewhitenoise; •, estimated levelshift(withitsestimated standard deviation in theparentheses directly underneath it); to,estimated trend;b'co, estimatedstandarddeviationof the estimatedtrend;O'&• ^* estimated"idealized" standard deviation of the estimated trend; r, the fraction of data before the intervention.
questionsabout changesto the environmentto be answered Since an approximate95% interval for log n* is given by more quickly.Efficient monitoringmay be used to show that logn* (•b) _+2•, an approximate 95%interval forn* is trends, which may have been expected, have not occurred. given by(•*e -2x/p,g*e2X/v).TheB givenin (4) is anestiSavingsfrom efficient monitoringwill include reducedmoni- mateof 2•/• by replacing the unknown 0 by •b.(Strictly toring costsaswell as reductionin environmentalimpact due speaking,this replacementintroducesadditional uncertainty. to early detection of changes. More sophisticatedinterval estimate based on the so-called variance-stabilizingtransformationis available, but the formula is more involved.) Appendix A1.
Data
Sources
A3.
Variance and Bias Calculations Under Models (5)
The UV data are from the 14 stations from the original and (1) Robertson-BergerUV monitoring network reported on by With Y = (Y•., -.., Yr)' as the T X 1 vector of observaWeatherhead et al. [1997]. The length of the UV data records tions, the trend model (5) with interventionlevel shift term rangefrom 10 to 18 years.The columnozonedata are version 7 from the total ozonemappingspectrometeron the Nimbus7 included may be expressedin matrix form as satellite.
The estimates were based on data from
The 14 locations
1979 to 1993.
for total ozone were chosen to coincide
with
Y:
+
the 14 UV monitoring stations.The surface ozone data are where X is a T x 3 matrix consistingof the constant,trend, from a network collected at the NOAA Climate Monitoring andinterventionterms,• = (/•, to,8)' representsthe vectorof and DiagnosticLaboratory;somestationsare analyzedby Olt- unknown regressioncoefficientsto be estimated, and N = mansand Levy [1994]. The length of the data recordsrange (N•,.", NT)' is the T x 1 vector of noiseterms. For the from 4 to 24 years.Data can be obtainedfrom S. Oltmans.The model above with AR(1) noise process {Nt}, let e = humidity, pressure, and solar radiation data are from the (V'i - 02N•.,e2, '", er) ', whichhasCov(e): o-•I.From National Renewable Energy Laboratory (NREL) data set the AR(1) equationN t - c•Nt_•. = et, it follows that the which includessomemodeled data not expectedto affect the noisevectorN satisfiesP'N = e, sothat N = P'-•-e, where the resultspresentedhere [NREL, 1992].The data are from years
matrixP' is T x T with(1, 1)-element equalto V'i - 02,the
1961 to 1991 inclusive. A2. Derivation in Section 2.4
for the Interval
Procedure
remaining diagonal elements equal to 1, the (i, i - 1)elementsequalto -0, and zero elementsotherwise.Hence the
Given
To derive a 95% confidenceinterval for n * when qbis esti-
covariance matrixof N hastheformCov(N) =Cov (P'-•e) = o-•P'-•-P -•- = o-•V,withV = P'-•-P-•- or V-•- = PP'. Consider
mated by•b,wefirstwriteexpli•citly n* = n*(0). Conse-the transformedequation quen•tly, whenqbis estimated by 0, n* is estimated byg* =
Y*: P'Y = P'X[5 + P'N = X*[5 + •, (A2) n*(0). To constructa 95% interval for n*, we will use a normal approximationon the log n* scale. The log scale is whereX* = P'X, and Cov(e) -- o-•I.The generalized least usedbecauseit helpsto improvethe normal approximationas squares(GLS) estimatorof [5 in model (A1) is then the ordiwell as to maintain the positivityof the resulting(interval)
naryleastsquares estimator in thetransformed model,• = estimate. Bythewell-known 8 method• (i.e.,first-order Taylor (X*'X*)-•-X*'Y *, withcovariance matrix expansion),the variancel/of log n *(0) is givenby
Cov(•) = c%2(X*'X*) -•.
•/-•-
xVar(•)=3(1-0)
M '
After some algebra,we have
(A3)
17,160
WEATHERHEAD
X*'X*
ET AL.: TREND DETECTION
h2 h3 hs ,
useful, simple approximationfor the variance expressionin (A5) is Var (&) • r% 2123/{(1 - qb)2r(r2 - 1)•, whichis
h4 hs h 6
quiteaccurate when
hih2 h4]
--
where
h 1= (T-
1
by assumingmodel (1), without interventionterm, when the true modelis (5) with a nonzerolevel shift interventionterm.
1)(1 - q•)2q- (1 - q•2),
Undermodel(5),E[Y*] = X•[• + U* 8, sotheexpected value
of • is
1--01
h2 •
12 [TT(T- 1)(1- qb) + T + qb],
=
1 1,U*8= [•l+ [hi h2 h hs 8. +(X,,X:)_•X, h:]-•[h4]
+ T2q•(]- q•)+ Tqb- q•2], h4: (T(T-
=
1
1
h3-- 1-•[• T(T + 1)(2T+ 1)(1- •)2
h5=
is not closeto 1.
Next,weconsider thebiasoftheGLSestimator • obtained
Hence we see from this that the expectedvalue of the trend estimateobtainedfrom model (1), when the true modelis (5)
T0)(1 - q•)2+ (1 _ q•),
with a nonzero level shift, is
T0)(1 -
[(T + T0)(1 - qb)+ 1 +
24
hihs-
+ &[To(To12 h 6 = (T-
[ro]= +
'
s.
In the caseof no autocorrelationwith qb= 0, this expression
T0)(1 - q•)2q_1.
reduces
to
Therefore usingmatrixalgebra, Var (&) isobtained aso• times the (2, 2)-elementof the inverseof the matrix X*'X*, which yields the expression 2
h2h4
h6(hlh6 -- h42)
6(T0-
E[&]= •o+ Now
Var(&)= 0-•(h•h6 - h42)(h3h6hs 2)- (h2h6h4h5) 2 _=rr•2h,2(qb, T, To).
(A4)
1)(T-
To + 1)
6r(1 - r)
r(r2 1)/12 8= •o+
we consider
the
situation
of section
3.2 where
8. it is
assumedthat for 8 in model(5) we haveadditionalinformation
in theformof a (prior)estimate 80withspecified variance •. Set [•o = (0, 0, 8o)'with 5;•- denotingthe 3 x 3 matrixwhose elementsare zero exceptfor the (3, 3)-elementwhichis equal
In the specialcaseof no autocorrelationwhen qb= 0, expres- to 1/•. Thentheefficientestimate of [• in (A1), whichincorsion (A4) simplifiesto poratesthe aboveprior estimateor additionalinformation,is
givenby• = (X*'X* + •E•-)-•(X*'Y * + o•E•-[•o) withco-
Var (&)
variance
matrix
123
2
•/3-- COV (•
= o-•2T(r+l)(2T+l)_3(r_ro+l)(r+ro)2_3ro2(ro_l) . A close approximation to the variance for this special case is Var(&) • rr•2123 / { r 3[1 - 3r(1 - r)]} =
o-•2/{n311 - 3r(1 - r)]}, wheren = T/12 is thenumberof years of data, and r = (T O - 1)/T is the fractionof preintervention data. In a similar way, for a general value of qb,a
usefulapproximationfor (A4) can be obtainedas Var (&)
o-•2/{ (1 - qb) 2n3[ 1 - 3r(1 - r) ] }, butthisapproximation is not so accuratefor larger values of
The resultfor model (1) can be obtainedas a specialcase, since(1) resultsfrom omittingthe level shift term BUt from (5). Let X = [X•, U] where X• representsthe T x 2 matrixof
-
[•) =o-•(X 2 ,, X ,
+
O.2•-) - 1.
(A6)
Notation
Y, St Xt Nt et /• •oo &
time seriesof environmentaldata. seasonalmean function. linear trend function. AR(1) noise. white noise. mean of data. trend to be detected. trend
estimate.
o-•, uncertaintyof trend estimate.
regressors formodel(1),and[5= ([5i, 8)' where[5•= (• •o)'. o-N standarddeviationof the AR(1) noise. Then the GLS estimator of [5• under model (1) is [5• = ..•_ , with covariancematrix
= O.e(XlX•)-i • ø'•h2 2 1 h
,
estimateof rrN. standarddeviationof the white noise. estimate of autocorrelationof the AR(1) noise.
qb estimate of
whereX• = P'X•. So,for model(1) it followsimmediately that Var (&) isgivenbyo• timesthe(2, 2)-element of theinverse of the matrixX•'X•, hi
2 • o.•2h 2( Var(&)= o-• h•h3 - h22 ok, T).
&N o-• a-• qb
(AS)
For the specialcaseof no autocorrelationwhen qb= 0, the expressionsimplifiesto Var (&) = r% 2123/{r(r2 - 1)). A
TO T n n*
time of interventionor level shift. length of data set in months. length of data set in years. yearsto detect a trend.
g*
estimate
of n*.
Acknowledgments. This researchwas supportedby the U.S. Departmentof Energyunder grant DE-FG03-94ER61857. We are grate-
WEATHERHEAD
ET AL.: TREND DETECTION
ful for their support.The assistance of Gene Maxwell in providingdata and adviceis also appreciated.
17,161
J. J. DeLuisi, C. L. Mateer, and D. J. Wuebbles, Effects of autocor-
relation and temporal samplingschemeson estimatesof trend and spatial correlation,J. Geophys.Res., 95, 20,507-20,517, 1990. Weatherhead, E. C., G. C. Tiao, G. C. Reinsel, J. E. Frederick, J. J.
References
DeLuisi, D. Choi, and W.-K. Tam, Analysisof long-termbehaviorof ultraviolet radiation measuredby Robertson-Bergermeters at 14 sitesin the United States,J. Geophys.Res., 102, 8737-8754, 1997. Zerefos, C. S., D: S. Balis, A. F. Bais, D. Gillotay, P. C. Simon, B. 70-79, 1975. Mayer, and G. Seckmeyer,Variability of UV-B at four stations in Box, G. E. P., G. M. Jenkins,and G. C. Reinsel, Time SeriesAnalysis: Europe, Geophys.Res.Lett., 24, 1363-1366, 1997. Forecastingand Control, 3rd ed., Prentice-Hall, Englewood Cliffs, Box, G. E. P., and G. C. Tiao, Intervention analysiswith applications to economicand environmentalproblems,J. Am. Stat. Assoc., 70,
N.J., 1994.
Krzygcin,J. W., UV controlling factors and trends derived from the ground-basedmeasurementstaken at Belsk, Poland, 1976-1994, J. Geophys.Res., 101, 16,797-16,805, 1996. Miller, A. J., R. M. Nagatani, G. C. Tiao, X. F. Niu, G. C. Reinsel, D. Wuebbles,and K. Grant, Comparisonsof observedozone and temperature trends in the lower stratosphere,Geophys.Res. Lett., 19, 929-932, 1992.
National RenewableEnergy Laboratory(NREL), National SolarRadiation Data Base (1961-1990) User'sManual, version 1.0, Golden, Colo., 1992.
Oltmans, S. J., and H. Levy II, Surface ozone measurementsfrom a global network,Atmos.Environ.,28, 9-22, 1994. Reinsel, G. C., G. C. Tiao, D. J. Wuebbles, J. B. Kerr, A. J. Miller,
R. M. Nagatani, L. Bishop,and L. H. Ying, Seasonaltrend analysis of published ground-basedand TOMS total ozone data through 1991,J. Geophys.Res., 99, 5449-5464, 1994. Scotto,J., G. Cotton, F. Urbach, D. Berger, and T. Fears, Biologically effective
ultraviolet
radiation:
Surface measurements
in the United
States, 1974 to 1985, Science,239, 762-764, 1988. Stolarski, R. S., P. Bloomfield, and R. D. McPeters, Total ozone trends
deducedfrom Nimbus 7 TOMS data, Geophys.Res.Lett., 18, 1015-
W.-K. Cheangand G. C. Reinsel,Department of Statistics,University of Wisconsin,Madison, WI 53706. D. Choi and X.-L. Meng, Department of Statistics,University of Chicago,Chicago,IL 60637. J. DeLuisi, NOAA ARL, 325 Broadway,Boulder, CO 80303. J. E. Frederick,Department of GeophysicalSciences,Universityof Chicago,Chicago,IL 60637. T. Keller, National Center for AtmosphericResearch,Boulder, CO 80307.
J. B. Kerr, AtmosphericEnvironment Service,4905 Dufferin Street, Downsview, Ontario, Canada.
A. J. Miller, NOAA NWS, 5200 Auth Road, Washington,D.C. 20233.
S. J. Oltmans,NOAA CMDL, 325 Broadway,Boulder, CO 80303. G. C. Tiao, Graduate School of Business,University of Chicago, Chicago,IL 60637. E. C. Weatherhead, CooperativeInstitute for Researchin the EnvironmentalSciences,Universityof Colorado,325 Broadway,Boulder, CO 80309. (e-mail:
[email protected]) D. J. Wuebbles,Departmentof AtmosphericSciences,Universityof Illinois, Champaign,IL 61801.
1018, 1991.
Stolarski,R., R. Bojkov, L. Bishop, C. Zerefos, J. Staehelin, and J. Zawodny, Measured trends in stratosphericozone, Science,256, 342-349, 1992. Tiao, G. C., G. C. Reinsel, D. Xu, J. H. Pedrick, X. Zhu, A. J. Miller,
(ReceivedNovember 10, 1997; revisedMarch 2, 1998; acceptedMarch 4, 1998.)
