Transcript Methodology - Meteorological Synthesizing Centre

Trend analysis: methodology
Victor Shatalov
Meteorological Synthesizing Centre East
TFMM trend analysis workshop, 17-18 November 2014
Main topics

Trend analysis of annual averages of
concentration/deposition fluxes

Trend analysis of monthly averages (with seasonal
variations)
TFMM trend analysis workshop, 17-18 November 2014
Trend analysis: generalities
Aim: investigation of general tendencies in time series such as:
 Measured and calculated pollutant concentrations at monitoring sites
 Average concentrations/deposition fluxes in EMEP countries
…
Method: trend analysis – decomposition of the considered series into
regular component (trend) and random component (residue)
Residue (random component)
1.2
1.0
1.0
0.8
0.8
Residue
Trend
0.4
ng/m3
0.6
0.4
0.2
0.2
0.0
0.0
-0.2
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
0.6
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
ng/m3
B[a]P concentrations in Germany
1.2
TFMM trend analysis workshop, 17-18 November 2014
Main steps
 Detection of trend and its character:
 increasing
 decreasing
 mixed
 Identification of trend type:
 linear
 quadratic
 exponential
 other
 Quantification of trend:
 total reduction
 annual reduction
 magnitude of seasonal variations
 magnitude of random component
 other
 Interpretation of the obtained results
Presentation by Markus Wallasch,
15 TFMM meeting, April 2014
TFMM trend analysis workshop, 17-18 November 2014
Determination of trend existence
B[a]P measurements: SE12
0.14
Decreasing
pair
0.12
0.10
0.08
0.06
Increasing
pair
0.04
0.02
B[a]P concentrations in Germany
Mixed trend character:
3
0.6
0.4
0.2
2010
2008
2006
2004
2002
2000
1998
1996
1994
0
1992
statistically significant (at 90% level)
increasing trend
Z = 1.8
0.8
1990
Typical
situation
for HMs
and –POPs
In the period
from 2004
to 2010
Z = - 4.05
1
ng/m
In the period from 1990 to 2000 –
statistically significant (at 95% level)
decreasing trend
1.2
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
Z = - 1.49
1997
0.00
1995
Decreasing trend at 85% significance
level
1996
ng/m3
Mann-Kendall test:
Z = (number of increasing pairs) –
(number of decreasing pairs) with
normalization.
Critical values:
± 1.44 at 85% level
± 1.65 at 90% level
± 1.96 at 95% level
TFMM trend analysis workshop, 17-18 November 2014
Determination of trend type: linear trend
Conc = A · Time + B + ω
ω – residues (random component)
Calculation of A and B:
regression or Sen’s slope
B[a]P concentrations in Germany
Random component
1.4
0.3
Calculations
1.2
Z = 3.8
increase
0.1
0.8
ng/m3
0.6
0.0
Residual trend exists
2010
2008
2006
2004
2002
2000
1998
1996
1994
2010
2008
2006
2004
2002
2000
1998
1996
-0.3
1994
0
1992
-0.2
1990
0.2
1992
-0.1
0.4
1990
3
1
ng/m
Z = - 3.1
decrease
0.2
Linear trend
Criterion of the choice of trend type: Mann-Kendall test should not show
statistically significant trend on all sub-periods of the time series
TFMM trend analysis workshop, 17-18 November 2014
Criterion of non-linearity
Criterion of non-linearity of the obtained trend in time:
chord)] · 100%
NL = max[abs(Δ
i
i /Ci
B[a]P concentrations in Finland
B[a]P concentrations in Germany
1.4
0.8
0.6
3
C
Trend
Δi
Chord
0.4
NL = 15.6%
Non-linear trend
Air concentrations
Trend
Fraction of non-linear trends
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
ng/m3
1.0
chord
i
ng/m
1.2
Air concentrations
0.100
0.090
0.080
0.070
0.060
0.050
0.040
0.030
0.020
0.010
0.000
0.2
Heavy metals
(Pb)
87%
B[a]P concentrations in Belgium
0.0
POPs
(B[a]P)
0.600
62%
NL = 8.1%
0.500
ng/m
3
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
0.700
0.400
0.300
0.200
0.100
Linear trend
Air concentrations
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
Supposed threshold value: 10%
Trend
0.000
TFMM trend analysis workshop, 17-18 November 2014
Determination of trend type: mono-exponential trend
Conc = A · exp(- Time / t) + ω,
Calculation of A and t:
least square method
t – characteristic time
Residual (random component)
B[a]P concentrations in Germany
0.3
1.4
Calculations
1.2
Exponential trend
ng/m3
0.6
0.1
0.0
0.4
-0.1
0.2
Residual trend exists
2010
2008
2006
2004
2002
2000
1998
1996
1994
1992
2010
2008
2006
2004
2002
2000
1998
1996
1994
1992
1990
-0.2
0
1990
ng/m
3
1
0.8
Z = 3.2
increase
Z = - 3.3
decrease
0.2
TFMM trend analysis workshop, 17-18 November 2014
Determination of trend type: polynomial trend
Conc = A ·
Time2
+ B · Time + C + ω
Calculation of A, B and C:
least square method
B[a]P concentrations in Germany
Random component
1.4
0.2
Calculations
1.2
Z = 0.5
no trend
Polynomial trend
0.1
ng/m3
0.8
0.6
0.4
0.0
-0.1
0.2
Residual trend exists
2010
2008
2006
2004
2002
2000
1998
1996
1994
1992
-0.2
1990
2010
2008
2006
2004
2002
2000
1998
1996
1994
1992
0
1990
3
1
ng/m
Z = -2.3
decrease
TFMM trend analysis workshop, 17-18 November 2014
Determination of trend type: bi-exponential trend
Conc = A1 · exp(- Time /t1) + A2 · exp(- Time /t2)
Ai – amplitudes, ti – characteristic times
Calculated by
least square method
Residual (random component)
B[a]P concentrations in Germany
0.2
1.4
Z=0
no trend
Calculations
1.2
Bi-exponential trend
0.1
ng/m3
0.8
0.6
0.4
0.0
-0.1
See [Smith, 2002]
2008
2006
2004
2002
2000
1998
1996
1994
statistically significant
residual trend obtained
1992
-0.2No
1990
2010
2008
2006
2004
2002
2000
1998
1996
1994
1992
0
2010
0.2
1990
3
1
ng/m
Z = -1.4
no trend
TFMM trend analysis workshop, 17-18 November 2014
Statistical significance of increasing trend
B[a]P concentrations in Germany
Typical situation for B[a]P:
increase in the end of the period
1.4
Calculations
1.2
Mann-Kendall test for 2004 – 2010:
3
does not confirm statistically significant
increasing trend
does not claim the absence of increasing
trend
0.8
0.6
0.4
0.2
1.2
1
1
0.8
2010
2008
2006
2004
2002
2000
1998
1996
1.2
0.8
3
3
ng/m
ng/m
0.6
0.4
0.6
0.4
0.2
2008
2006
2004
2002
2000
1998
1996
1994
1992
Increase
is statistically
0
significant
1990
2010
2006
2004
2002
2000
1998
1996
1994
1992
-0.2
2010
0.2
0
1990
[A, B] – confidence
interval for slope of
random component
1994
1992
1990
0
Confidence interval for trend slope:
[TS0 + A, TS0 + B]
TS0 – slope of
calculated trend
Z = 1.8
2008

1
ng/m

Bi-exponential trend
TFMM trend analysis workshop, 17-18 November 2014
Non-linear trend analysis
Conc = A1 · exp(- Year / t1) + A2 · exp(- Year / t2) + ω
Regression model, non-linear in the parameters t1 and t2
Non-linear regression models are widely investigated, for example:
 Nonlinear regression, Gordon K. Smith, in Encyclopedia on
Environmetrics, ISBN 0471899976, Wiley&Sons, 2002, vol 3, pp.
1405 – 1411
 Estimating and Validating Nonlinear Regression Metamodels in
Simulation, I. R. dos Santos and A. M. O. Porta Nova,
Communications in Statistics, Simulation and Computation, 2007,
vol. 36: pp. 123 – 137
 Nonlinear regression, G. A. F. Seber and C. J. Wild, WileyInterscience, 2003
TFMM trend analysis workshop, 17-18 November 2014
Parameters for trend characterization:
reduction/growth
B[a]P concentrations in Germany
1.4
Calculations
1.2
ng/m3
Cbeg
Bi-exponential trend
Total reduction per period
Rtot = (Сbeg–Cend)/Cbeg=1–Cend/Cbeg
1.0
0.8
ΔCi
0.6
0.4
Relative annual reduction
Ri = ΔCi / Ci = (1 – Ci+1 / Ci)
Cend Average annual reduction
0.2
Rav = 1 – (Cend / Cbeg) 1/(N-1)
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
0.0
where N – number of years
Negative values of reduction mean growth
Reduction parameters
Rmin = min (Ri)
Rmax = max (Ri)
Rav
Rtot
For the considered example:
Rmin = - 6% (growth)
Rmax = 15%
Rav = 6%
Rtot = 69%
TFMM trend analysis workshop, 17-18 November 2014
Parameters for trend characterization:
random component
Normalized
random
component
Residue
(random
component)
B[a]P concentrations in Germany
30% 0.4
1.4
Calculations
1.2
Bi-exponential trend
20%
0.3
1.0
0.4
Δ
0.1
0%
0.0
0.2
Frand
1990
1991
1990
1992
1991
1993
1992
1994
1993
1994
1995
1995
1996
1996
1997
1997
1998
1998
1999
1999
2000
2000
2001
2001
2002
2002
2003
2003
2004
2004
2005
2005
2006
2006
2007
2007
2008
2008
2009
2009
2010
2010
0.6
ng/m3
10%
-10%
-0.1
0.0
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
ng/m3
0.2
0.8
-20% -0.2
Parameter: standard deviation of random component
normalized by trend values Frand = σ(Δ/Ctrend)
For the considered example: Frand = 11%
TFMM trend analysis workshop, 17-18 November 2014
Seasonal variations of pollution
B[a]P concentrations at CZ3 site
3.5
3.0
2.5
ng/m3
B[a]P concentrations measured
at EMEP site CZ3 from 1996 to
2010. Pronounced seasonal
variations are seen.
2.0
1.5
1.0
0.5
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
0.0
Pb concentrations at DE7 site
60
50
ng/m3
Pb concentrations measured at
EMEP site DE7 from 1990 to
2008. Seasonal variations are
also seen.
70
40
30
20
10
Seasonal variations are characteristic of heavy metals and (particularly)
for POPs
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
1990
0
TFMM trend analysis workshop, 17-18 November 2014
Possible approaches to description of seasonal variations
t – time
Bi-exponential approximation
Conc = A1 · exp(– t / t1) · (1 + B1 · cos(2p · t – φ1)) +
A2 · exp(– t / t2) · (1 + B2 · cos(2p · t – φ2))
t – chatracteristic times,
A, B – constants,
φ – phase shifts.
Mono-exponential approximation *)
Conc = A · exp(– t / t + B · cos(2p · t – φ))
or
Log(Conc) = A’ – t / t + B · cos(2p · t – φ)
*) Kong et al., Statistical analysis of long-term monitoring data…
Environ. Sci. Techn., 10/2014
TFMM trend analysis workshop, 17-18 November 2014
Usage of higher harmonics
Measurement data at CZ3 from 1996 to 2010
3.5
Measurements
3
trend standard
Trend calculated by
bi-exponential approach.
Possible artifact:
negative trend values
ng/m
3
2.5
2
1.5
1
0.5
0
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
-0.5
Possibility to avoid negative values: usage of higher harmonics
Conc = Tr1 + Tr2 ,
Tri = Ai·exp(– t / ti)·(1+Bi·cos(2p·t–φi)+Ci·cos(4p·t–ψi))
3.5
Measurements
3
two harmonics
2.5
1.5
1
0.5
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
0
1996
ng/m
3
Statistical significance
of second harmonic:
Fisher’s test F
2
TFMM trend analysis workshop, 17-18 November 2014
Usage of higher harmonics
Average B[a]P concentrations in Europe from 1990 to 2010
(main harmonic only)
One harmonic
0.35
Calculations
0.3
Trend
Main component
0.25
ng/m
3
Poor approximation
for small values of
concentrations
0.2
0.15
0.1
0.05
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
1990
0
Residues for one-harmonic approximation
Residues, one harmonic
0.1
0.08
Pronounced
harmonic trend with
doubled frequency
0.04
0.02
0
-0.02
-0.04
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
-0.06
1990
ng/m
3
0.06
TFMM trend analysis workshop, 17-18 November 2014
Usage of higher harmonics
Average B[a]P concentrations in Europe from 1990 to 2010
(main harmonic only)
One harmonic
0.35
Calculations
0.3
Trend
Main component
0.25
ng/m
3
Poor approximation
for small values of
concentrations
0.2
0.15
0.1
0.05
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
1990
0
Trend including two harmonics
0.35
Calculations
Trend
Main component
0.3
Significance of
second harmonic is
confirmed by
Fisher’s test
0.2
0.15
0.1
0.05
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
0
1990
ng/m
3
0.25
TFMM trend analysis workshop, 17-18 November 2014
Splitting trends to particular components
Full trend
2
Example: average B[a]P
concentrations for Germany
from 1990 to 2010.
Concentrations
1.8
Trend
Main component
Ctot
1.6
1.4
ng/m
3
1.2
1
0.8
0.6
0.4
Ctot = Cmain + Cseas + Crand
0.2
2009
2010
2008
2009
2008
2006
2004
2007
2005
2003
2004
2010
2007
2006
2005
2002
2001
2000
1998
1997
1996
1995
1994
1993
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
-0.6
1995
0
1994
-0.4
1993
0.2
1992
-0.2
1992
0
0.4
1990
2003
0.2
1991
0.6
1990
3
ng/m
0.8
2002
2001
2000
1999
1998
1997
Crand
0.4
1999
Cmain
1991
3
1996
0.8
0.6
1
ng/m
1995
Random component
Main component
1.2
Full trend
1994
1993
1992
1991
1990
0
Seasonal component
0.8
Cseas
0.6
Relative annual reductions
(as above): Rmin, Rmax, Rav, Rtot
0.4
0
-0.2
-0.4
-0.6
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
-0.8
1990
ng/m
3
0.2
TFMM trend analysis workshop, 17-18 November 2014
Splitting trends to particular components
Full trend
2
Example: average B[a]P
concentrations for Germany
from 1990 to 2010.
Concentrations
1.8
Trend
Main component
Ctot
1.6
1.4
ng/m
3
1.2
1
0.8
0.6
0.4
Ctot = Cmain + Cseas + Crand
0.2
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
Fraction of0.6trends with essential
seasonality
C
rand
0.4
93%
3
Heavy0.20metals (Pb)
ng/m
0.6
0.4
-0.2
POPs-0.4(B[a]P)
0.2
100%
Seasonal
component
Seasonal
component,
normalized
Cseas
80%
0.6
60%
0.4
40%
0.2
20%
0%0
Average value of the annual amplitude
of the normalized seasonal
component Fseas
ng/m
3
0.8
100%
-20%
-0.2
-40%
-0.4
-60%
-0.6
-80%
Normalization: Cseas/Cmain
2010
2010
2009
2009
2008
2008
2007
2007
2006
2006
2005
2005
2004
2004
2003
2003
2002
2002
2001
2001
2000
2000
1999
1999
1998
1998
1997
1997
1996
1996
1995
1995
1994
1994
1993
1993
1992
1992
1991
1991
1990
1990
-0.8
-100%
Threshold value: 10%
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
1990
1991
-0.6
0
1990
3
ng/m
Full trend
1996
0.8
Cmain
0.8
1995
Random component
Main component
1.2
1
1994
1993
1992
1991
1990
0
TFMM trend analysis workshop, 17-18 November 2014
Splitting trends to particular components
Full trend
2
Example: average B[a]P
concentrations for Germany
from 1990 to 2010.
Concentrations
1.8
Trend
Main component
Ctot
1.6
1.4
ng/m
3
1.2
1
0.8
0.6
0.4
Ctot = Cmain + Cseas + Crand
0.2
Cseas
0.6
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1997
1996
1998
2010
2010
2009
2009
2008
2008
2007
2007
2006
2006
2005
2005
2004
2004
2003
2003
2002
2002
2001
2001
2000
2000
1999
1999
1998
1998
1997
1997
1996
1996
1995
1995
1994
1994
1993
1993
Normalization: Crand/Cmain
Seasonal component
0.8
1992
1992
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
1990
0
1991
1991
0.2
1990
1990
0.6
0.4
0.4
Standard deviation of normalized
random component Frand
3
0.2
0
-0.2
-0.4
-0.6
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
-0.8
1990
ng/m
1995
3
20%
0.2
0%
0
-20%
-0.2
-40%
-60%
-0.4
-80%
-0.6
-100%
Crand
ng/m
3
ng/m
Full trend
1994
0.8
100%
80%
0.6
60%
0.4
40%
Cmain
0.8
1993
Random
component,
normalized
Random
component
Main component
1.2
1
1992
1991
1990
0
TFMM trend analysis workshop, 17-18 November 2014
Phase shift as a fingerprint of source type
14
Anthropogenic
Secondary
Trends for PB concentrations at CZ1
Air concentrations, ng/m
3
12
10
8
6
4
Δφ
2
2010
2009
2008
2007
2006
2005
2004
2003
2002
2001
2000
1999
1998
1997
1996
1995
1994
1993
1992
1991
1990
0
Difference Δφ of phase shift φ between Pb pollution at CZ1 location due to
anthropogenic and secondary sources.
Phase shift can be used to determine which source type (anthropogenic or
secondary) mainly contributes to the pollution at given location (in a particular
country).
TFMM trend analysis workshop, 17-18 November 2014
List of trend parameters
Parameters for trend characterization:
 Relative reduction over the whole period (Rtot),
 Relative annual reductions of contamination:
 average over the period (Rav),
 maximum (Rmax),
 minimum (Rmin).
 Relative contribution of seasonal variability (Fseas).
 Relative contribution of random component (Frand).
 Phase shift of maximum values of contamination with respect to the
beginning of the year (φ).
Statistical tests:
 Non-linearity parameter (NL)
 Relative contribution of seasonal variability (Fseas)
10%
10%