FORECASTING CONDITIONAL CLIMATECHANGE USING A HYBRID
APPROACH
By
AKBAR AKBARI ESFAHANI
B.S., University of Colorado Denver, 2010
A thesis submitted to the
Faculty of the Graduate School of the
University of Colorado in partial fulfillment
of the requirements for the degree of
Masters of Integrated Sciences
Integrated Sciences
2013
This thesis for the Masters of Integerated Science degree by
Akbar Akbari Esfahani
has been approved for the
Integrated Sciences Program
By
Weldon Lodwick, Chair
Loren Cobb
Michael J. Friedel
September 23, 2013
Akbari Esfahani, Akbar (MSIS, Integrated Sciences)
Modeling ClimateChange in the Southwestern United States
Thesis directed by Professor Weldon Lodwick.
ABSTRACT
A novel approach is proposed to forecast the likelihood of climatechange across spatial
landscape gradients. This hybrid approach involves reconstructing past precipitation and
temperature using the selforganizing map technique; determining quantile trends in the
climatechange variables by quantile regression modeling; and computing conditional
forecasts of climatechange variables based on selfsimilarity in quantile trends using the
fractionally differenced autoregressive integrated moving average technique. The
proposed modeling approach is applied to states (Arizona, California, Colorado, Nevada,
New Mexico, and Utah) in the southwestern US, where conditional forecasts of climate
change variables are evaluated against recent (2012) observations, evaluated at a future
time period (2030), and evaluated as future trends (20092059). These results have broad
economic, political, and social implications because they quantify uncertainty in climate
change forecasts affecting various sectors of society. Another benefit of the proposed
hybrid approach is that it can be extended to any spatiotemporal scale providing self
similarity exists.
The form and content of this abstract are approved. I recommend its publication.
Approved: Weldon Lodwick
DEDICATION
I dedicate this work to the tireless support of my loving wife Naimeh and to my
parents who made the biggest sacrifice by leaving behind their family so that we could
have a better life.
IV
ACKNOWLEDGMENTS
I would like to thank Michael J. Friedel, my teacher, mentor,
without his tireless push, this research would have not
and friend. For
been possible
v
TABLE OF CONTENTS
Chapter
Table...........................................................................vii
Figure.........................................................................viii
1. Introduction..................................................................1
2. Methodology...................................................................4
2.1 Reconstructing climatechange variables...............................4
2.2 Quantile trends in climatechange variables...........................7
2.3 Forecasts of climatechange variables.................................8
2.3.1 Test of determinism..............................................8
2.3.2 Determining the fractal dimension...............................10
2.3.3 Longmemory process.............................................11
2.3.4 Fractionally differenced autoregressive integrated moving average
model.................................................................12
2.4 Software for reconstruction, quantile trends, and forecasting and empirical
cumulative distribution functions.........................................13
3. Results......................................................................14
3.1 Reconstructing climatechange variables..............................14
3.2 Quantile trends in climatechange variables..........................20
3.3 Forecasts of climatechange variables................................21
3.3.1 A look at the past..............................................21
3.3.2 A look at the present...........................................23
3.3.2 A look at the future............................................30
4. Conclusions..................................................................46
References.......................................................................47
vi
TABLES
Table
3.1  Summary of longmemory Hurst parameters for precipitation and temperature in
Arizona, California, Colorado, New Mexico, Nevada, and Utah...............................22
3.2  Forecast trends (20002050) for temperature, precipitation, and drought in Arizona,
California, Colorado, New Mexico, Nevada, and Utah. Results are for 0.25, 0.50, and
0.75 quantiles............................................................................43
Vll
FIGURES
Figure
2.1 Schematic depicting the hybridmodeling framework used to forecast climate
change................................................................................4
3.1  Quantile (upper panel) modeling of Palmer Drought Severity Index (PDSI)
measurements (02009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect
0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The
approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue) are
indicated by double arrows............................................................16
3.2  Quantile (upper panel) modeling of reconstructed temperature measurements (0
2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25
(yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the
Medieval Warming Period (red) and Little Ice Age (blue) are indicated by double arrows.
17
3.3  Quantile (upper panel) modeling of reconstructed precipitation measurements (0
2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25
(yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the
Medieval Warming Period (red) and Little Ice Age (blue) are indicated by arrows.......19
3.4  Comparison of conditional forecasts and measurement observations for 2012 Palmer
Drought Severity index (PDSI): (a) California, (b) Colorado, and (c) Nevada. Arrows
denote median forecasted (thin) and observed (open) PDSI value........................26
3.5  Comparison of conditional forecasts and measurement observations for 2012
temperature: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median
forecasted (thin) and observed (open) temperature value...............................28
3.6  Comparison of conditional forecasts and measurement observations for 2012
precipitation: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median
forecasted (thin) and observed (open) precipitation value.............................29
3.7  Conditional forecasts and histogram (lower panel) for Palmer Drought Severity
Index in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote
median forecasted 2030 (thin) and observed 2012 (open) Palmer Drought Severity Index
values................................................................................32
3.8  Conditional forecasts and histogram (lower panel) for temperature in the year 2030:
(a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030
(thin) and observed 2012 (open) temperature values....................................36
3.9  Conditional forecasts and histogram (lower panel) for precipitation in the year 2030:
(a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030
(thin) and observed 2012 (open) precipitation values..................................40
viii
1. Introduction
People benefit from a multitude of resources and processes supplied by natural
ecosystems (Randhir and Ekness, 2009). These benefits include water resources suitable
for supporting various sectors of society such as agriculture, construction, daily living,
energy, fishing, forestry, manufacturing, public health, recreation, transportation, and
overall economic development that maintains life systems prompting sustainability.
Climate change is frequently cited as one external driver of ecosystems (Furnis, 2010).
Because climate is temporally and spatially dependent, change at a global scale differs
from regional or local scales (Friedel, 2012b). One reason for spatial differences is the
superposition of largescale climate patterns due to atmospheric and oceanic
teleconnections (Schwing et al., 2002). Climate change also differs across temporal scales
over which there are variations in amplitude, gradient, and duration (Woodhouse and
Overpeck, 1998; Mann et al., 2008; Mann et al., 2009; Friedel, 2011a, Friedel, 2012a).
In many studies, the duration of climate change is considered shortterm (years to
decades) variability (Woodhouse and Overpeck, 1998, Ubilava and Helmers, 2013). Short
term climate variability is often attributed to oscillations in the sea surface temperature
(SST) that alter ocean currents and overlying air pressure resulting in a redistribution of
temperature and precipitation (Smith and Reynolds, 2003). The El Nino Southern
Oscillation (ENSO) is considered the strongest shortterm periodic fluctuation (27 years)
with a rise (El Nino) or decrease (La Nina) of SST in the equatorial Pacific Ocean (Blade
et al., 2008, Ubilava and Helmers, 2013). The influence of this teleconnection is not
uniform in the United States, and ENSO events can affect things like water supply, water
quality, riparian habitat, power generation, and range productivity. Related drought
1
consequences often include crop failure, debris flows, insect infestations, pestilence,
violent conflict, wildfires, and disruptions to economic and social activities (Riebsame et
al., 1991). Longterm climate variability (hundreds to thousands of years) is often
attributed to alterations in geologic and extraterrestrial processes, such as volcanic
aerosols (Rampino and Self, 1982) and solar activity (Gray et al., 2010). Longterm
climatechange reconstructions provide insight on past surface temperature and drought
variability at timescales crossing centuries or millennia (Mann et al., 2009; Friedel,
2012a). Should natural or anthropogenic forcing influence the frequency or intensity of
climate change, there is an increased likelihood for future ENSO hazards placing national
and global security at risk (Riebsame et al., 1991).
For these reasons, climatechange forecasts could benefit many sectors of society,
but the scaledependent complexities render it a challenge using traditional processbased
models. Specifically, climate forcing is known to interact with ecosystems characterized
by coupled, nonlinear, and multivariate processes. Data associated with these ecosystems
2
are typically sparsely populated ranging spatially from local (1000s km ) to global and
temporally from immediate (l10s years) to longterm (100s to 1000s years). One critical
issue is the lack of essential calibration data that results in large inaccuracies (Loke et al.,
1999). Other nonprocessbased modeling efforts include linear timeseries models, such
as the Autoregressive (AR) or Autoregressive Integrated Moving Average (ARIMA)
models (Said and Dickey, 1984; Cowpertwait and Metcalfe, 2009). These traditional
linear timeseries modeling schemes are too rigid with respect to detecting unexpected
features like the onset of nonlinear trends, or patterns restricted to subsamples of a data
set.
2
One alternative paradigm is to use a hybrid approach in which soft computing
provides data for subsequent use in traditional numerical or empirical modeling. Some
applications of the hybrid modeling approach are in rainfallrunoff (Jain and Kumar,
2007), debrisflow (Friedel, 2011b), mineralresource (Friedel, 2012b), and unexploded
ordinance (Friedel et. al., 2012). In this study, the goal was to evaluate the efficacy of
hybrid modeling for forecasting climatechange over spatial landscape gradients in the
southwestern United States. The objectives were to: (1) evaluate independent
observations against the 2012 forecast of temperature, precipitation, and drought across
California, Colorado, and Nevada; (2) evaluate future forecast trends over 50 years of
probable temperature, precipitation, and drought across Arizona, California, Colorado,
New Mexico, Nevada, and Utah; and (3) to evaluate the performance characteristics of
models that are generated (Bennett et al, 2013). However, in todays modeling world, it is
not just enough to evaluate the performance of a model, but it is a necessity to account for
model uncertainty. It is to this end that this study introduces a new modeling paradigm to
account for uncertainty in timeseries forecasting.
This study extends the work of Friedel (2011a) who sought to reconstruct 2,000
years of past temperature and precipitation for the southcentral and southwestern United
States; Esfahani and Friedel (2010) who identified a longmemory process in
reconstructed climate variables; Friedel (2012a) who used quantile regression to quantify
uncertainty in global reconstructions of past temperature and precipitation and Caballero
et al. (2002) and Nunes et al. (2011) who used a fractal approach to investigate the long
memory process associated with temperature.
3
2. Methodology
In modeling climatechange variables, the hybrid approach relies on four
computational steps: reconstruction, trends, forecasts, and uncertainty (Fig. 2.1). Each of
these steps is briefly described in the following sections.
Fig. 2.1 Schematic depicting the hybridmodeling framework used to forecast climate
change.
2.1 Reconstructing climatechange variables
The reconstruction of climatechange variables (temperature and precipitation)
follows the approach described by Friedel (2012a). In that approach, a selforganizing
map (SOM) technique was used to project input data to a discrete lattice of competitive
4
neurons (Kohonen, 2001), and estimates past climatechange variables by minimizing
topological error vectors (Fraser and Dickson, 2007). The process of projecting data is
essentially a data compression technique (Hastie et al., 2002) for which the success of
topologypreservation was analyzed based on the quantization error E(G, X), given by:
M
E(G,X) =^2j2jhijjxj Wi2, (1)
ieQ j=i
where w* are weight vectors assigned to a fixed number of N neurons in the map grid G,
Xj are the M input data vectors (economic mineralresource variables), hUI is a
neighborhood function, xy Wj is the Euclidian norm, and / is the best matching unit
(BMU) vector.
Implementation of the SOM learning method is based on the stochastic gradient
described by Kohonen (2001). It consists of a twostep process that is performed each
time an input pattern is presented to the map: competition to determine the BMU and
cooperative learning (spreading information contained in the current input vector across
the map). At the beginning of the unsupervised training phase, the weight vectors are
initialized to small random numbers. The input data vectors are presented to the map grid
in a random fashion to generate data clusters without introducing bias for a specific class.
In the first step, the BMU with map coordinates (/,, Ij) is determined as the grid neuron,
whose weight vector is the closest to the input given by:
I = ArgMini jeGx(i) w(j), (2)
where ArgMin is the minimum distance defining the central position of the neighborhood
function. The neighborhood function h, i is chosen to be Gaussian function given by:
5
(3)
hu(n) = exp
cr(n);
where r; ri corresponds to the distance between map neuron r, and BMU in the map
grid, and O(n) defines the width of the neighborhood function, a monotonically
decreasing function of the iteration (also called epoch) number n. In the second step, a
weight update is determined which is a function of the distance to the current BMU, as
expressed through the neighborhood function h,i(n). The weights are gradually adjusted
according to:
Wj(n + 1) = Wj(n) + a(n)hu(n)[x(n) W;(n)], (4)
where a(n) is a scalar value called the learning rate bounded on the interval [0, 1], The
BMU ensures that the largest weight correction (hiti(n)=1) is adjusted in the direction of
the input vector. The association effect takes place at the neighboring nodes but to a
lesser degree because of the Gaussian shape. This adaption procedure stretches the
weight vectors of the BMU and its topological neighbors towards the input vector.
Presenting similar input vectors to the map provides further activations in the same
neighborhood and thereby tends to produce clustering of data in the feature space.
Association between neurons decreases during the learning process (the width of the
neighborhood function &(n) is forced to decrease with n preserving large clusters of data
while enabling the separation of clusters that are closely spaced). Ultimately, this training
process results in a topology where similarities among data patterns are mapped into
similar weights of the neighboring neurons, and the asymptotic local density of the
weights approach that of the training set (Ritter and Schulten, 1986).
6
Crossvalidation (Efron and Tibshirani, 1993) is conducted to ensure the SOM
provides unbiased estimates of climatechange variables. In this case, known data values
are estimated based on distances among the available model vectors (Wang, 2003; Kalteh
and Berndtsson 2008). In the traditional approach, estimates of values are taken directly
from the prototype vectors of the best matching units (Fessant and Midenet, 2002; Wang,
2003). Often times certain training data sets result in biased estimates (Dickson and
Giblin, 2007; Malek, et al., 2008) requiring a modified scheme that incorporates
bootstrapping (Breiman, 1996), ensemble average (Rallo et al., 2002), or nearest neighbor
(Malek et al., 2008). This study uses an alternative iterative estimation scheme that
minimizes the topological error vector (Fessant and Midenet, 2002). The estimation of
past climatechange values for all variables is done simultaneously and referred to here as
the reconstruction. For more details about SOM training and estimation, the reader also is
referred to (Kohonen, 2001; Vesanto, and Alhoniemi, 2000).
2.2 Quantile trends in climatechange variables
The determination of quantile trends in the reconstructed climatechange variables
follows the approach described by Friedel (2012). One advantage of this approach is its
flexibility in modeling data with conditional functions that may have systematic
differences in dispersion, tail behavior, and other covariate features (Koenker, 2005). In
adopting this approach, the quantile regression approach is adopted thereby reducing the
ddimensional nonparameteric regression problems to a series of additive univariate
problems given by:
7
n n2
min
aeM
^pi(y; i) A^Jd^a
(5)
i=l i=l
where min is the minimization operator and at estimated coefficients at distinct xi. dj is
an ndimensional vector with elements (hjl, (hj^ + hj1), hj^) in the j, j+1, and j+2
positions and zeros elsewhere. ht = xi+l xt; pt is the quantile regression loss function;
A is a factor that controls the degree of smoothing; T is the transpose; xi are response
observation values associated with random variable X; and yÂ£ are response observation
values associated with random variable Y.
Individual quantile curves can then be specified as a linear bspline of the form:
where Yt_1) is the rth conditional quantile function, /?Â£(r) are the regression
quantiles; /?Â£ = ai+1 aijhi, are regression coefficients; 0Â£(.):i = 1 denote the
the spline have been selected, such models are then linear in parameters and therefore can
be estimated (Koenker, 2005).
2.3 Forecasts of climatechange variables
2.3.1 Test of determinism
Before fractal techniques can be used to analyze the reconstructed climatechange
data, the determinism of a fractal time series must be established. The data must be
identified as deterministic or statistical (Kaplan, 1994; Turcotte, 1997). The other
assumption to be satisfied is that stationarity exists for longmemory process data (Beran,
p
(6)
i=l
basis function of the spline; and r is the quantile (0 < r < 1). Once the knot positions of
8
1994). The determinism test developed by Kaplan (1994), allows the determination as to
whether a time series has a deterministic or stochastic structure. Turcotte (1997) defines
the set as a deterministic fractal set, if the set is scale invariant at all scales; and a
statistical fractal set, if the set is different at different scales but the differences do not
allow the scale to be determined. This is an important distinction because it determines
the correct way to calculate the fractal dimensions (as defined by equation (10)) for a
time series.
Mathematically, the deltaepsilon test of continuity (Kaplan, 1994) is applied to
determine the deterministic structure of a time series. To do so, the continuity test is
applied to orbits comprising the phase space topology created by timedelayed
embedding of the original set. This process is facilitated by generating an ensemble set of
surrogate time series using a bootstrapping approach (Constantine et al., 2010). Next, the
phasespace statistic, called the Estatistic (Kaplan, 1994), is calculated for the time
delayed embedding of the original time series as well as the ensemble of surrogate data.
The structure is then judged by comparing the Estatistics of the original set to the E
statistic of the surrogate set. A separation between the two statistics implies the existence
of a deterministic structure, and the converse implies the set is a realization of a random
process and thus has a stochastic structure.
According to Kaplan (1994), the EStatistic is defined as:
(7)
(8)
(9)
Sj,k = z(j) z(k),
Â£j,k = Kj + k) z(k+K),
s(r) = for j, k s. t. r < 5j k < r + Ar,
9
where 5j k is the Euclidean distance between phase space points z(j) and z(k), and ej k is
the corresponding separation distance between the points at a time k points in the future
along their respective orbits. The variable k is the orbital lag and the future points are the
images of the original pair. The increment Ar is the width of a specified Euclidean bin
size. Given Ar, the distance 5j k is used to identify the proper bin in which to store the
image distance ej k and the average of each bin forms the s(r) statistic. Finally, the E
statistic is formed by calculating the cumulative sum of the s(r) statistic (Constantine and
Percival, 2011).
2.3.2 Determining the fractal dimension
Fractal sets are characterized by their dimensions (Mandelbrot, 1967; Mandelbrot
and van Ness, 1968). Fractal dimensions are mathematically different from topological
dimensions in that while topological dimensions are strictly constructed from integers,
fractal dimensions can be fractional. While there are an infinite number of fractal
dimensions, the fractal information dimension is important because it describes the
entropy of a data set. In general, the fractal dimension is described by:
1 dlog^p?
q 1 d log(r)
(10)
where q refers to the dimension of interest. By embedding the dataset in an ndimensional
grid with cells having sides of size r, we compute the frequency for which a data point
falls into the ith cell, p off The information dimension is described by lim^i Dq. As the
information dimension approaches 1, the numerator of equation (4) changes to the
Shannons entropy describing changes in the entropy and trend of a time series (Barbara,
10
1999). Under these conditions, the fractal dimension of a stochastic selfsimilar data set
simplifies to:
1 logLp?
Dq q 1 log(r)
(11)
Without the partial differential equation requirement the problem becomes
computationally tractable and can be calculated easily, while there are several ways to
calculate the fractal dimensions, for the purposes of this paper, the box counting
technique is used (De Pison et al., 2008).
2.3.3 Longmemory process
The existence of a longmemory process was first explored by Hurst (1951) when
trying to find a solution on how to regularize flow of the Nile River. He observed that
long periods of highflow levels were followed by long periods of lowflow levels
indicating the existence of a longmemory process. In general, longmemory process
describes the process where by a time series shows that past values have a strong effect
on present or future values. Mandelbrot and van Ness (1968) introduced the Hurst
parameter (H) to describe the longterm memory of a timeseries process. The Hurst
parameter is given by:
H=l (12)
where A is the longmemory parameter with the range of 0 < A < 1; thus, H has a range
of < H < 1. The closer H is to 1, the more persistent the time series is considered and
at values less than or equal to a longmemory process does not exist. Currently there
are several methods to estimate the Hurst parameter; however, Rae et al. (2011) showed
11
that Whittles maximum likelihood estimation approach produced the best results (Beran,
1994).
2.3.4 Fractionally differenced autoregressive integrated moving average model
If climatechange variables have a longmemory process then the fractionally
differenced autoregressive integrated moving average model (FARIMA) can be applied
and forecasts generated. The FARIMA model is similar to the BoxJenkins ARIMA
model (Cowpertwait and Metcalfe, 2009) except that the integrated part of the ARIMA
model can be a fractional number defined (Beran, 1994) as:
d = H (13)
This expression underscores the need for a good estimation of the Hurst parameter
when dealing with climate timeseries data. If the climatechange data exhibits long
memory process and stationarity, the FARIMA model assumptions are satisfied and can
be used to generate forecasts following three steps. First, the model parameters for each
time series are determined using the maximumlikelihood estimators (Fraley et al., 2011).
Second, the model is fitted and forecasts made using the estimated parameters (Hyndman,
2011). Third, the forecasts are validated using various statistical measures (Hyndman and
Koehler, 2006).
2.3.5 Empirical cumulative distribution function
The use of quantile regression facilitates quantization of the prediction uncertainty
by constructing empirical cumulative distribution functions (ECDF) at any year of a
forecast. The empirical cumulative distribution function (ECDF) is determined by
modeling the collection of quantile forecast results. The ECDF is a step function fv
12
(14)
71
Fn( t) =
(Xj < t) 1
1 V"
n
i=i
where i is the number of tied observations at that value (missing values are ignored). For
observations x = (x1( x2, ..., xn), and Fn is the fraction of observations less or equal to t.
2.4 Software for reconstruction, quantile trends, and forecasting and empirical
cumulative distribution functions
The data mining, reconstruction, and analysis are carried out using the SiroSOM
(CSIRO Exploration & Mining, 2008) graphical user interface (GUI). This GUI provides
an interface between data sets and functions in the freely available SOM Toolbox
(Adaptive Informatics Research Center, 2010). Quantile modeling is conducted using the
quantreg, splines, and stats packages (Hornik, 2011); fractal transformation and fractal
dimension calculations are conducted using the fdim package (Martinez de Pison,
Ascacibar et al. 2008). Fractionally differenced autoregressive integrated moving
average modeling is conducted using the forecast and fractal package (Constantine and
Percival, 2011, Hyndman 2011) and the empirical cumulative distribution function
modeling is conducted using the stats package (Hornik, 2011); freely available at
http://www.rproject.org/ in the R toolbox.
13
3. Results
3.1 Reconstructing climatechange variables
Annual temperature and precipitation values were reconstructed across a gradient
of modern climate zones: Arizona (Desert), California (Mediterranean), Colorado
(Semiarid to Alpine), Nevada (Semiarid to Arid), New Mexico (Semiarid), and Utah
(Semiarid). The simultaneous reconstruction of annual temperature and precipitation was
done based on the selforganized nonlinear datavector relations among approximately
2,000 years (0 to 2009 AD) of reconstructed warmseason (average of June, July, and
August) Palmer Drought Severity Index (PDSI) data (Cook et al., 2004), 114 years
(18952009) of annual state precipitation (accumulation over January through
December) and temperature (average of January through December) data (National
Climatic Data Center, 2010), and other related tropical and extratropical measurements
described by Friedel (2012a). The temperature and precipitation are standard climate
variables, whereas the PDSI defines annual dry, neutral, and wet periods based on tree
ring information (Palmer, 1965; Cook et al., 2004). For a comprehensive review of the
reconstruction and validation of these data, the reader is referred to the cited references.
In this study, the reconstructions were verified against independent precipitation and
temperature data for the years: 1896, 1900, 1911, 1919, 1923, 1935, 1940, 1952, 1960,
1966, 1968 (La Nina), 1986, 1993, 1998 (La Nina), and 2005 (El Nino) using split and
crossvalidation (leave one out) approaches (Bennett et al, 2013). The Spearman Rho
correlation among observed and reconstructed values was greater than 95% with a p
value of 0.001. The magnitude of climate variability over the past 2,000 years is visually
apparent when inspecting plots for PDSI (Fig. 3.1 a, b, and c), temperature (Fig. 3.2 a, b,
14
PDSI PDSI
and c), and precipitation (Fig. 3.3 a, b, and c), where (a) is California, (b) is Colorado,
and (c) is Nevada.
15
0 500 1000 1500 2000
Year
(c) Nevada
Fig. 3.1 Quantile (upper panel) and histogram (lower panel) modeling of Palmer
Drought Severity Index (PDSI) measurements (02009): (a) California, (b) Colorado,
and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue),
0.95 (purple) quantiles. The approximate timing of the Medieval Warming Period (red)
and Little Ice Age (blue) are indicated by double arrows.
(a) California
Year
16
(b) Colorado
Year
(c) Nevada
Year
Fig. 3.2 Quantile (upper panel) and histogram (lower panel) modeling of reconstructed
temperature measurements (02009): (a) California, (b) Colorado, and (c) Nevada. The
trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles.
The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue)
are indicated by double arrows.
17
Precipitation, mm/year
250 300 350 400 450 500
(a) California
Precipitation, mm/year
400 600 800 1000
o
o
LO
i_ CO
CO
a)
c O
E CO
c
o
V'
ro
Q.
O
0)
O
LO
CM
O
O
CM
O
LO
Year
(c) Nevada
Fig. 3.3 Quantile (upper panel) and histogram (lower panel) modeling of reconstructed
precipitation measurements (02009): (a) California, (b) Colorado, and (c) Nevada. The
trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles.
The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue)
are indicated by arrows.
In general, the plots reveal that all states in the southwest experienced past
conditions ranging from extreme drought to extremely moist. It is interesting to note that
over the past 2000 years, the nodrought condition has the greatest frequency of
occurrence in California (62.5%), Colorado (58.5%), and Nevada (55%). This suggests
that while there is the likelihood for additional extreme events, there is a greater
likelihood for current drought conditions to be eventually mitigated. Regarding
temperature, the respective mode values in California, Colorado, and Nevada are about
14.8 C (range from 14 C to 16 C), 6.75 C (range from 5.5 C to 9.5 C), and 9.75 C
(range from 8 C to 12 C). Whereas the respective mode values for precipitation in
California, Colorado, and Nevada are about 425 mm (range from 150 mm to 1050 mm),
19
350 mm (range from 230 mm to 540 mm), and 190 mm (range from 120 mm to 350 mm).
Both the temperature and precipitation for these states reflect their association with the
modern climate gradient.
3.2 Quantile trends in climatechange variables
Quantile regression modeling is applied to the past (02009) PDSI (Fig. 3.1 a, b,
&c), temperature (Fig. 3.2 a, b, &c), and precipitation (Fig. 3.3 a, b, &c) data. The annual
quantiles trends are presented for climate change data in a) California, b) Colorado, and
c) Nevada. The decadal trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue),
0.95 (purple) quantiles determined using 200 bspline degrees of freedom. In general, the
quantile trends reveal that the longterm regional climate was interrupted by shortterm
changes. The socalled Medieval Warm Period (900 to 1250) and Little Ice Age
(1300 to 1850) appear as primary lower frequency disruptions over the last two
millennia in California and Nevada, and secondary (muted) disruptions in Colorado
(Crowley and Lowery, 2000; Mann, 2002). The approximate timing of the Medieval
Warming (horizontal red arrow) coincides with a decrease in PDSI (drier conditions),
increase in temperature, and decrease in precipitation; whereas the Little Ice Age Periods
(horizontal blue arrow) coincides with an increase in PDSI (wetter conditions), decrease
in temperature, and increase in precipitation. These findings are attributed to strong
ENSO teleconnections with California and Nevada, but mixed ENSO signals in Colorado
(it is a region between El Nino and La Nina latitudes). In computing forecasts for future
climate in southwestern states, additional quantile trends (0.05, 0.10, 0.20, 0.30, 0.40,
0.50, 0.60, 0.70, 0.80, 0.90, and 0.95) are used to increase resolution when modeling the
empirical cumulative distribution functions.
20
3.3 Forecasts of climatechange variables
3.3.1 A look at the past
To evaluate past climatechange data for longmemory process (Barbara, 1999),
the PDSI, temperature, and precipitation are transformed to an equivalent fractal
information dimension. The transformation is applied to 20year intervals resulting in
time series plots comprising 100 lags. The various methods used to calculate the long
memory Hurst parameter (H) are summarized in table 3.1. Based on these H calculations,
the precipitation appears to have the strongest longmemory process, whereas
temperature and PDSI have a comparatively weaker longmemory process. This finding
implies that the autocorrelation of precipitation is stronger than temperature; that is, the
precipitation persists over a longer number of years than temperature.
Another finding is that the arithmetic calculation of H falls within the bounds of
the Whittle approximation. One exception may be calculations based on the spectral
regression method. For example, application of this method to precipitation results in an
H parameter value that is outside the bounds of the Whittle approximation. This suggests
that the spectral regression method might not provide reasonable approximations for
comparatively strong longmemory processes. It is also interesting to note that the two
driest states, Nevada and Arizona, have the lowest H values for precipitation. This
suggests that wet cycles in these states are of shorter duration than other southwestern
states. By contrast, Colorado has the second shortest temperature cycle and California the
longest temperature cycle. At the same time, Colorado has the longest drought cycle and
California has the shortest drought cycle. These findings suggest that Colorado is more
likely to experience extended drought cycles at different temperatures, whereas
21
California is more likely to experience longer cycles of high temperatures that do not
result in droughts, which can be explained by the fact that California borders the Pacific
Ocean.
Table 3.1 Summary of longmemory Hurst parameters for precipitation and
temperature in Arizona, California, Colorado, New Mexico, Nevada, and Utah.
Precipitation
Ha calculation H interval Whittle Method H via Spectral Regression
2.50% 97.50%
AZ 0.948 0.850 1.171 0.693
CA 0.978 0.869 1.234 0.701
CO 0.979 0.871 1.190 0.742
NM 0.983 0.930 1.368 0.742
NV 0.964 0.832 1.152 0.693
UT 0.982 0.934 1.312 0.999
Temperature
Ha calculation H interval Whittle Method H via Spectral Regression
2.50% 97.50%
AZ 0.606 0.569 0.651 0.616
CA 0.662 0.546 0.790 0.597
CO 0.579 0.533 0.628 0.592
NM 0.602 0.521 0.738 0.600
NV 0.560 0.521 0.611 0.570
UT 0.605 0.647 0.895 0.601
Palmer Drought Severity Inc ex
Ha calculation H interval Whittle Method H via Spectral Regression
2.50% 97.50%
AZ 0.639 0.590 0.696 0.647
CA 0.554 0.507 0.604 0.538
CO 0.651 0.596 0.707 0.682
NM 0.633 0.574 0.692 0.649
NV 0.630 0.572 0.688 0.646
UT 0.598 0.533 0.664 0.638
3.3.2 A look at the present
The previous section established that the southwest climatechange
reconstructions are stationary and stochastically selfsimilar. Given these facts, the
22
FARIMA model is applied to individual quantile trends for forecasting of climatechange
variables. In all cases, the fitted FARIMA models are characterized as one of the
following types: from AR(0) MA(0) to AR(0) MA(6). Using the fitted FARIMA models,
forecasts are generated over the range of quantiles from which ECDFs are then computed
for each variable.
The FARIMA model performance is evaluated by comparing the median 2012
state forecasts of PDSI (Fig. 3.4), temperature (Fig. 3.5), and precipitation (Fig. 3.6)
(California, Colorado, and Nevada respectively) to state observations (National Oceanic
and Atmospheric Administration, 2013). In this comparison, the PDSI reflects values
averaged over the period of JuneJulyAugust, the temperature represents values
averaged over the period of JanuaryDecember, and the precipitation reflect values
accumulated precipitation over the period of JanuaryDecember. When comparing the
2012 climatechange observations to forecasts, it is important to note the differences in
their spatiotemporal representation. Specifically, the observed PDSI and reconstructed
climate variables are values associated with a grid location (influenced by a 250 km x
250 km region) defined by Cook et al. (2004), whereas the observed climate variables
represent values averaged over climate divisions with larger areas defined by the National
Oceanic and Atmospheric Administration (2013).
The respective 2012 state observations and forecasts of PDSI, temperature, and
precipitation are presented in Figs. 3.43.6. Reviewing the plots for California, Colorado,
and Nevada reveals that all observations plot within the probable forecast limits
supporting the usefulness of the hybrid approach. In addition, the range of probable
climate variable differs among states reflecting the modem climate gradient and
23
underscoring the stationarity of conditions over the past 2000 years. For example, the
respective forecast range of PDSI (Fig. 3.4) for California, Colorado, and Nevada is from
3.50 to 0.5, from 4.25 to 2.25, and from 5.0 to 1.0. The respective forecast range of
temperature (Fig. 3.5) for California, Colorado, and Nevada is from 15.1 C to 15.8 C,
from 6.65 C to 8.0 C, and from 9.6 C to 10.5 C. The respective forecast range of
precipitation (Fig. 3.6) for California, Colorado, and Nevada is from 200 mm to 625 mm,
from 340 mm to 440 mm, and from 170 mm to 260 mm.
Comparing the median forecast (indicated by vertical arrows associated with the
0.5 quantile) to observed climate variables averaged over different climate divisions
reveals reasonable correspondence. For example, the respective median forecast and
averaged observations of PDSI in California (Fig. 3.4a) are about 2.50 and 2.20
(climate divisions 4 and 6). The respective median forecast and averaged observation of
PDSI in Colorado (Fig. 3.4b) are about 0.75 and 2.71 (climate divisions 4 and 1); and
the respective median forecast and averaged observation of PDSI in Nevada (Fig. 3.4c)
are about 2.15 and 2.17 (climate divisions 3 and 4). Similarly, the respective median
forecast and averaged observation of temperature in California (Fig. 3.5a) are about 15.44
C and 15.2 C (climate divisions 4 and 6). The respective median forecast and averaged
observation of temperature in Colorado (Fig. 3.5b) are about 7.42 C and 7.78 C (climate
divisions 4 and 1); and the respective median forecast and averaged observation of
temperature in Nevada (Fig. 3.5c) are about 10.05 C and 10.6 C (climate divisions 3 and
4). The respective median forecast and averaged observation of precipitation in California
(Fig. 3.6a) are about 450 mm and 440 mm (climate divisions 4 and 6). The respective
median forecast and averaged observation of precipitation in Colorado (Fig. 3.6b) are
24
about 403 mm and 392 mm (climate divisions 4 and 1); and the respective median
forecast and averaged observation of precipitation in Nevada (Fig. 3.6c) are about 202
mm and 188.5 mm (climate divisions 3 and 4). Comparison of median forecasts to values
for individual divisions and all state divisions reveals the heterogeneous nature of climate
variables within each state; that is, the observations across climate regions in California
are mostly Mediterranean; in Colorado, they range from semiarid to alpine; in Nevada
they range from semiarid to arid. This finding suggests that future improvements,
potentially spanning five to 15 years, to forecasting may be achieved by introducing
PDSI data associated with grid nodes across all state climate divisions.
Severe Moderate No
Drought Drought Drought
3.00 2.00 1.99
to to to
3.99 2.99 + 1.99
4 1 1 1 r Ibservation climate Observation climate division 4 all divisions (1.66) (2.52)
3 2 1 0
(a) California
PDSI
25
Extreme
Drought
<4
Severe
Drought
3.00
to
3.99
Moderate
Drought
2.00
to
2.99
No
Drought
1.99
to
+1.99
Moderately
Moist
+2.00
to
+3.99
PDSI
(b) Colorado
Extreme Severe Moderate No
Drought Drought Drought Drought
<4 3.00 2.00 1.99
to to to
3.99 2.99 +1.99
PDSI
(c) Nevada
Fig. 3.4 Comparison of conditional forecasts and measurement observations for 2012
Palmer Drought Severity index (PDSI): (a) California, (b) Colorado, and (c) Nevada.
Arrows denote median forecasted (thin) and observed (open) PDSI value.
26
quantiles quantiles
O
CO
o
(Â£>
O
O
CM
O
O
O
Observation climate
divisions 4 and 6
(15.2 C)
&
...........i i
i i i i r
15.1 15.2 15.3 15.4 15.5
Observation climate
divisions all
(15.1 C)
i r~
15.6 15.7
Temperature, Celsius
(a) California
Observation climate
divisions 4
(7.78)
Temperature, Celsius
(b) Colorado
27
Observation climate
divisions 3
(10.6)
9.6
9.8
10.0 10.2
Temperature, Celsius
Observation climate
divisions all
(9.89 C)
10.4
(c) Nevada
Fig. 3.5 Comparison of conditional forecasts and measurement observations for 2012
temperature: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median
forecasted (thin) and observed (open) temperature value.
Observation climate 1
divisions 4 and 6
(440 mm)
r Observation climate
r divisions all
rJ (494 mm)
r r
300 400 500 600
Precipitation, mm/year
(a) California
28
(/>
0)
c
TO
Z5
O"
o
CO
o
CD
O
o
CN
O
O
O
Observation climate
divisions 1 and 4
(392 mm)
Observation climate
divisions all
(385.8 mm)
360 380 400 420
Precipitation, mm/year
(b) Colorado
tr
180
200 220
Precipitation, mm/year
240
260
Observation climate
divisions 3 and 4
(147 mm)
(c) Nevada
Fig. 3.6 Comparison of conditional forecasts and measurement observations for 2012
precipitation: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median
forecasted (thin) and observed (open) precipitation value.
29
3.3.2 A look at the future
The conditional state 2030 forecasts are presented together with 2012
observations (average of climate divisions) for PDSI, temperature, and precipitation
(Figs. 3.73.9). For example, the conditional forecasts of PDSI in California, Colorado,
and Nevada are presented in Fig. 3.7.
5 0 5
Extreme Moderate No Moderately Extremely
Drought Drought Drought Moist Moist
<4 2.00 1.99 +2.00 >+4
to to to
2.99 + 1.99 + 2.99
(a) California
30
Frequency quantiles
0 5
Moderate No Moderately Extremely
Drought Drought Moist Moist
2.00 1.99 +2.00 >4
to to to
2.99 +199 +3.99
(b) Colorado
31
5 0 5
Extreme Moderate No Moderately Extremely
Drought Drought Drought Moist Moist
<4 2.00 1.99 +2.00 >+4
to to to
2.99 +1.99 +3.99
(c) Nevada
Fig. 3.7 Conditional forecasts and histogram (lower panel) for Palmer Drought
Severity Index in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows
denote median forecasted 2030 (thin) and observed 2012 (open) Palmer Drought
Severity Index values.
In California (Fig. 3.7a), the median PDSI forecast value of about 1.42 is larger
than the observed 2012 value of about 2.20 (average of climate divisions 4 and 6). This
32
suggests the likelihood for a shift toward moister conditions; that is, the current moderate
drought condition is likely to shift to a midrange condition within the next 8 years. In
Colorado (Fig. 3.7b), the median PDSI forecast value of about 0.40 is larger than the
observed 2012 value (average of climate divisions 4 and 1) of about 3.50, also
suggesting a shift from severe drought to a midrange condition in the next 8 years. In
Nevada (Fig. 3.7c), the median PDSI forecast value of about 1.60 is larger than the
observed 2012 value of about 2.17 (average of climate divisions 4 and 6), also
suggesting the likelihood for a shift toward moister conditions in the next 8 years. In
addition to trends based on median values, there also is some probability that California
could experience conditions ranging from severe drought to moderately moist in the next
8 years. Colorado could experience conditions ranging from extreme drought to
moderately moist in the next 8 years and Nevada could experience conditions ranging
from severe drought to midrange in the next 8 years. Relative to Colorado, California
tends toward the moist and Nevada shifts toward drought conditions.
33
Frequency quantiles
14.0 14.5 15.0 15.5 16.0
Temp in Bins of 0.5 degrees of Celcius
(a) California
34
Frequency quantiles
2012 Observation
(7.78)
5 6 7 8 9
(b) Colorado
Temp in Bins of 0.5 degrees of Celcius
2012 Observation
(10.6)
Temp in Bins of 0.5 degrees of Celcius
(c) Nevada
Fig. 3.8 Conditional forecasts and histogram (lower panel) for temperature in the year
2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted
2030 (thin) and observed 2012 (open) temperature values.
The conditional state 2030 forecasts of temperature in California, Colorado, and
Nevada are presented in Fig. 3.8. In California (Fig. 3.8a), the median temperature
36
forecast value of about 15.28 C is about the same as the observed 2012 temperature
(average of climate divisions 4 and 6) of about 15.2 C. This suggests the likelihood for a
slight shift toward warmer conditions in the next few years. In Colorado (Fig. 3.8b), the
median temperature forecast value of about 7.3 C is less than the 2012 temperature
(average of climate divisions 4) of about 7.78 C, suggesting a shift toward cooler
conditions in the next few years. In Nevada (Fig. 3.8c), the median temperature forecast
value of about 9.9 C is smaller than the 2012 temperature (average of climate divisions 4
and 6) of about 10.6 C, suggesting the likelihood for cooler conditions in the next few
years. In addition to trends based on median values, there also is some probability that
California could experience temperatures ranging from about 14.7 C to 15.8 C,
Colorado could experience temperatures ranging from about 6.5 C to 8.25 C, and
Nevada could experience temperatures ranging from about 9.2 C to 10.7 C. Within the
next two to eight years, relative to Colorado, California and Nevada tend toward warmer
conditions.
37
Frequency quantiles
200 400 600 800 1000
Precip in Bins of 50 mm/year
(a) California
38
Frequency quantiles
250 300 350 400 450 500 550
Precip in Bins of 50 mm/year
(b) Colorado
39
100 150 200 250 300 350 400
Precip in Bins of 50 mm/year
(a) Nevada
Fig. 3.9 Conditional forecasts and histogram (lower panel) for precipitation in the year
2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted
2030 (thin) and observed 2012 (open) precipitation values.
The conditional state 2030 forecasts of precipitation in California, Colorado, and
Nevada are presented in Fig. 3.9. In California (Fig. 3.9a), the median precipitation
40
forecast value of about 490 mm is larger than the 2012 precipitation value (average of
climate divisions 4 and 6) of about 440 mm. This suggests the likelihood for a future shift
toward drier conditions in the next few years. In Colorado (Fig. 3.9b), the median
precipitation forecast value of about 390 mm is about the same as the 2012 precipitation
value (average of climate divisions 4) of about 392 mm, suggesting no change in
precipitation conditions. In Nevada (Fig. 3.9c), the median temperature forecast value of
about 205 mm is larger than the historical median value of about 188.5 mm, suggesting
the likelihood for future wetter conditions in roughly eight years. In addition to trends
based on median values, there also is some probability that California could experience
precipitation amounts ranging from 275 mm to 725 mm, Colorado could experience
precipitation amounts ranging from 300 mm to 470 mm, and Nevada could experience
precipitation amounts ranging from 150 mm to 300 mm within the next eight years.
Relative to Colorado, California tends toward wetter conditions and Nevada toward drier
conditions.
In addition to evaluating forecasts for the year 2030, an analysis of trends is
conducted over forecasts spanning a 50year period (from 2010 to 2060). A summary of
findings is presented for the southwestern states in Table 2. These results indicate that the
southwestern US is likely to experience decreasing temperatures that will likely continue
over the next 25 years resulting in wetter and cooler conditions for the region. The decay
represents a decrease of 0.20 degrees Celsius over 20 years in temperature. The accuracy
of the forecast results for temperatures are verified by comparing the root mean squared
error (RMSE) (Bennett et al, 2013; Hyndman and Koehler, 2006) of the forecast of each
41
state to the mean of the forecast of that state. For temperature, the range is from 2.3% to
8.2% of RMSE to the mean indicating a low variability.
The PDSI index will increase in the same period, indicating wetter conditions
across the southwestern states. However, this increase is only by about 0.8 on the PDSI
index, which suggests only a very moderate decrease of drought. The exception to this
increase is California where the forecast suggests a constant level of PDSI suggesting that
California has reached its equilibrium, the end of a longmemory process. Since PDSI
values are mainly negative, RMSE is not a good indicator of performance, since it relies
on the absolute values, instead the mean error (ME) of the residuals is chosen here. A
comparison between the ME of each state and the mean of the forecast of each state
indicates variability of about 2% to 4%.
42
Table 3.2 Forecast trends (20002050) for temperature, precipitation, and drought in
Arizona, California, Colorado, New Mexico, Nevada, and Utah. Results are for 0.25,
0.50, and 0.75 quantiles.
Palmer Drought Severity Index (PDSI)
Forecast Accuracy Mean of Series Forecast Direction
State ME1
AZ 0.020 0.54 increase*
CA 0.016 0.24 constant
CO 0.011 0.51 . increase
NM 0.014 0.61 . increase
NV 0.018 0.70 . increase
UT 0.009 0.45 . increase
Increase = increase of 0.8 to 1.0 on PDSI scale over Decrease = decrease of 0.8 to 1.0 on PDSI scale over Since PDSI is measured from positive to negative, increase indie 0 years 10 years ates wetter conditions
Precipitation
Forecast Accuracy Mean of Series Forecast Direction
State RMSE2 MAE3
AZ 59.5 44.2 306.1 sharp increase
CA 98.3 77.0 528.2 increase
CO 45.7 35.7 387.7 decrease
NM 55.6 41.9 333.8 increase
NV 31.3 24.9 215.0 decrease
UT 39.4 31.1 283.9 increase
Increase = increase of 5 mm over 10 years Sharp Increase = increase of 20 mm over 10 years Decrease = decrease of 5 mm over 10 years Sharp Decrease = decrease of 20 mm over 10 years
Temperature
Median of Series Forecast Accuracy Mean of Series Forecast Direction
RMSE MAE
AZ 0.62 0.51 15.6 decrease
CA 0.39 0.32 15.1 decrease
CO 0.59 0.47 7.2 decrease
NM 0.46 0.37 11.8 decrease
NV 0.58 0.45 9.7 decrease
UT 0.60 0.50 9.0 decrease
Decrease = decrease of 0.20 degrees over 20 years Increase = increase of 0.20 degrees over 20 years
1 ME mean error, since PDSI is mostly a negative value, MAE and RMSE do not apply here.
2 RMSE root mean squared error.
3 MAE mean absolute error.
43
The forecast results for precipitation is slightly different than those of temperature
and PDSI. The forecast accuracy is measured using the RMSE for each state, which
indicates a variability of 11% to 19%. The higher variability can be explained by the fact
that precipitation has much stronger longmemory cycles (150 to 300 years) and a fifty
year forecast is only part of the cycle. It is interesting to note that Nevada and Colorado
will experience an increase in precipitation in the near future (next 5 to 10 years),
however; the long term forecast horizon (10 to 20 years) calls for less precipitation. This
indicates that the short forecast horizon for these 2 states is at the end of one cycle and
the beginning of the next cycle.
44
4. Conclusions
The proposed hybrid modeling approach is useful for forecasting drought,
temperature, and precipitation at the level of state climate divisions for a span of 1 to 15
years. This process requires the reconstruction of past climate variables that are stationary
and stochastically selfsimilar with quantile regression modeling used to facilitate
quantization of forecast uncertainty. The application to southwestern states provided
reconstructions of past (02009) climate records exhibiting temporal paleoclimatic
features, such as the Medieval Warm Period and Little Ice Age. Independent performance
testing using modem (2012) observations averaged over appropriate climate divisions
demonstrated good correspondence to median forecasts.
This finding supports the 50 years of forecasting (20102060) to assist managers
in formulating decisions, potential mitigating strategies, and policy associated with
future, uncertain climatechange. Differences among median forecasts and observations
from other climate divisions demonstrate the heterogeneous nature of climate variables
within each state. This finding supports future improvements in forecasting by
introducing additional paleoclimatic data associated with grid nodes crossing all of the
state climate divisions. Because the proposed hybrid approach can be extended to any
spatiotemporal scale providing selfsimilarity exists, forecasting has the possibility to
address economic, political, and social aspects affecting various sectors of society. As
with any datadriven approach, the introduction of additional related information can be
expected to further reduce the uncertainty.
While the approach presented was geared toward climate change, the work here
can be extended to any type of time series analysis. It represents a superior method to
45
traditional accuracy measures of forecasting. That is, traditional measures are tuned to
account for errors in the model, however; the method presented here accounts for the
uncertainty of the entire modelling process. The method can be applied to time series and
spatial analysis in fields such as economics, ecology, and biology. Beyond the novel
method presented for forecasting uncertainty, the fractal information dimension
transformation of data can be used to find changes to time series. Current research in the
area concetrates on the correlation of volcanic acitivities to climate shocks found in the
fractal transformation of 2000 years of data presented in this research (Friedel and Akbari
Esfahani, 2013).
46
REFERENCES
Akbari Esfahani, A., Friedel, M. J. 2010. The fractal nature of climate change 2000
years in retrospect, H21G1126 presented at 2010 Fall Meeting, AGU, San
Francisco, Calif., 1317 Dec.
Akbari Esfahani, A., Friedel, M.J., 2012, A Fractional Forecast to the Climate of
Southwestern United States, Mathematics of Planet Earth presented at 2012
Meeting, IUGG Conference on Mathematical Geophysics, Edinburgh, United
Kingdom, 1822 June.
Adaptive Informatics Research Center (AIRC), 2010. SOM Toolbox, Helsinki University
of Technology, Laboratory of Computer and Information Science, Adaptive
Informatics Research Center, http://www.cis.hut.fi/projects/somtoolbox/
Barbara, D., 1999. Chaotic mining: Knowledge discovery using the fractal dimension,
George Mason University Information and Software Engineering Department,
Fairfax, VA.
Bennet, N.D., Croke, F.W., Guariso, G., Guillaume, J.H.A., Hamilton, S.H., Jakeman,
A.J., MarsiliLibelli, S., Newham, L.T.H., Norton, J.P., Perrin, C., Pierce, S.A.,
Robson, B., Seffelt, R., Voinov, A.A., Fath, B.D., Andressian, V., 2013,
Characterising performance of environmental models, Environmental Modelling &
Software 40 (2013), 120.
Beran, J., 1994. Statistics for longmemory process, Boca Ranton: CRC Press LLC.
Blade, I., Newman, M., Alexander, M.A., Scott, J.D., 2008. The late fall extratropical
response to ENSO: Sensitivity to coupling and convection in the tropical west
Pacific. J. Climate 21(23), 61016118.
47
Breiman, L., 1996. Bagging predictors, Machine Learning, 24, 123140.
Caballero, R., Jewson, S. Brix, A., 2002. Long memory in surface air temperature:
Detection, modeling, and application to weather derivative valuation. Climate
Research 21, 127140.
Constantine W., Percival, D., 2010. fractal: Insightful Fractal Time Series Modeling and
Analysis. R package version 1.03. http://CRAN.Rproject.org/package=fractal
Cook, E. R., Woodhouse, C. A. Meko, D. M. Stahle, D. W., 2004. Long term aridity
changes in the western United States, Science 306 (5698), 10151018.
Constantine, W., Percival, D., 2011. Fractal: fractal time series modeling and analysis. R
package version 1.11, http://CRAN.Rproject.org/package=fractal, 2011.
Cowpertwait, P.S., Metcalfe, A.V., 2009. Introductory time series with R, New York:
Springer.
Crowley, T.J., Lowery, T.S., 2000. How Warm Was the Medieval Warm Period?
AMBIO: A Journal of the Human Environment 29 (1), 5154.
CSIRO Exploration & Mining, 2008, SIROSOM software.
http://www.csiro.au/science/SOMdataanalysis.html
Dickson, B.L., Giblin, A.M., 2007. An evaluation of methods for imputation of missing
trace element data in groundwaters. Geochemistry: Exploration, Environment, and
Analysis 7 (2), 173178.
Efron B., Tibshirani R.J., 1993. An introduction to the bootstrap, in Monographs on
Statistics and Applied Probability, 57, 436 pp, Chapman & Hall, London.
Fessant, F., Midenet, S., 2002. Selforganizing map for data imputation and correction in
surveys. Neural Computing Applications 10, 300310.
48
Fraley, C. Leisch, F., Maechler, M., Reisen, V., Lemonte, A., 2011. Fracdiff:
Fractionally differenced ARIMA aka ARFIMA(p,d,q) models. http://CRAN.R
proj ect. org/package=fracdiff
Fraser S., Dickson B.L., 2007. A new method for data integration and integrated data
interpretation: selforganizing maps. In Proceedings of Exploration 07: Fifth
Decennial International Conference on Mineral Exploration, edited by B. Milkereit,
907910.
Friedel, M.J., 2011a. Climate change effects on ecosystem services in the United States
Issues of national and global security. In: Baba, A., Tayfur, G., Howard, K.W.F.,
Friedel, M.J., Chambel, A., 2011, Climate Change and its Effect on Water
Resources Issues of National and Global Security, NATO Science for Peace and
Security Science Series C. Environmental Security, vol. 3, Springer, Dordrecht,
The Netherlands, 318 p. ISBN 9789400711426.
Friedel, M.J., 2011b. A datadriven approach for modeling postfire debrisflow volumes
and their uncertainty, Environmental Modelling & Software 26 (12), 15831598.
Friedel, M.J., 2012a, Datadriven modeling of surface temperature anomaly and solar
activity trends, Environmental Modelling & Software 37, 217232.
Friedel, M.J., 2012b, Hybrid modeling to predict the economic feasibility of mining
undiscovered porphyry copper deposits. Applied Soft Computing 13, 10161032.
Friedel., M.J., Asch, T., Oden, C., 2012. Hybrid analysis of multiaxis electromagnetic
data for discrimination of munitions and explosives of concern, Geophysical
Journal International 190 (2), 960980.
49
Friedel, M. J. Akbari Esfahani, A., 2013. Influence of global volcanic processes on
climate modulation in the southwestern US, T24 presented at 2013 Fall Meeting,
GSA, Denver, Colorado., 28 Oct. 1 Nov.
Furnis, M.J., 2010. Water, climate change, and forests. U.S. Forest service, General
Technical Report, PNWGTR812, pp. 75.
Gray, L.J., Beer, J., Geller, M., Haigh, J.D., Lockwood, M., Matthes, K. Cubasch, U.,
Fleitmann, D., Harrison, G., Hood, L., Luterbacher, J., Meehl, G.A., Shindell, D.,
van Geel, B., White, W., 2010. Solar influences on climate. Reviews of Geophysics
48: RG4001, 153.
Haslett, J., Raftery, A.E., (1989). Spacetime Modelling with Longmemory Dependence:
Assessing Ireland's Wind Power Resource (with Discussion); Applied Statistics 38,
150.
Hastie, T., Tibshirani, R., Friedman, J., 2002. The elements of statistical learning,
SpringerVerlag, Berlin, 533 pp.
Hornik, R., 2011. The R FAQ, ISBN 3900051089. http://www.rproject.org/, last
visited 1192012.
Hurst, H.E., 1951. Long term storage capacity of reservoirs, Trans. Am. Soc. Civ. Eng.,
116, 770779, 1951
Hyndman, R.J., 2011. Forecast: Forecasting functions for time series. http://CRAN.R
project.org/package=forecast.
Hyndman, R.J., Koehler, A.B., 2006. Another look at measures of forecast accuracy,
International Journal of Forecasting 22 (4), 679688.
50
Jain, A., Kumar, A.M., 2007. Hybrid neural network models for hydrologic time series
forecasting. Applied Soft Computing 7 (2), 585592.
Kalteh, A.M., Berndtsson, R., 2008. Review of the selforganizing map (SOM) approach
in water resources: Analysis, modeling and application. Environmental Modelling
& Software 23, 835845.
Kaplan, D.T., 1994. Exceptional events as evidence for determinism. Physica D 73, 38
48.
Koenker, R., 2005. Quantile regression, Econometric Society Monographs, Cambridge
University Press, New York, 349 p.
Kohonen T., 1984. SelfOrganization and Associative Memory. Springer, Berlin.
Kohonen, T., 2001. Selforganizing maps, Third Extended Edition, Springer Series in
Information Sciences, Vol. 30, Springer, Berlin, Heidelberg, New York, 253 p.
Loke, E., ArnbjergNielsen, K., Harremoes, P., 1999. Artificial neural networks and grey
box modelling: a comparison. In: Joliffe, I.B., Ball, J.E. (Eds.), Eighth International
Conference: Urban Storm Drainage Proceedings, vol. 1. The Institution of
Engineers Australia, Australia.
Malek, M.A., Harun, S., Shamsuddin, S.M., Mohamad I., 2008. Imputation of time series
data via Kohonen selforganizing maps in the presence of missing data. Engineering
and Technology 41, 501506.
Mandelbrot, B.B., 1967. How Long is the coast of Britain? statistical selfsimilarity and
fractional dimension, Science 156: 636638, 1967.
Mandelbrot, B.B., van Ness, J.W., 1968. Fractional brownian motions, fractional noises
and applications, Siam Review 10 (4), 422437.
51
Mann, M.E., 2002. Little Ice Age. In Volume 1, The Earth system: Physical and chemical
dimensions of global environmental change, Encyclopedia of Global Environmental
Change (pp. 504509). Chichester: John Wiley & Sons.
Mann, M. E., Zhang, Z., Hughes, M. K., Bradley, R. S., Miller, S. K., Rutherford, S., Ni,
F., 2008. Proxybased reconstructions of hemispheric and global surface
temperature variations over the past two millennia. Proceedings of the National
Academy of Sciences, 105 (36), 1325213257.
Mann, M. E., Zhang, Z., Rutherford, S., Bradley, R. S., Hughes, M. K., Shindell, D., Ni,
F., 2009. Global signatures and dynamical origins of the little ice age and medieval
climate anomaly. Science 326: 12561260, 2009.
Martinez de Pison Ascacibar, F.J., Ordieres Meres, Joaquin, Castejon Limas, Manuel, de
Cos Juez, F.J., 2008, Functions for calculating fractal dimensionPackage fdim,
http://cran.rproject.org/src/contrib/Archive/fdim/
National Oceanic and Atmospheric Administration, 2013. National Climatic Data Center.
http://vlb.ncdc.noaa.gov/tempandprecip/timeseries/
Nunes, S.A., Romani, L.A., Avila, A.M., Traina, C.J., de Sousa, E.P., Traina, A.J., 2011.
Fractalbased analysis to identify trend changes in multiple climate time series.
Journal of Information and Data Management 2 (1), 5157.
Palmer, W.S., 1965. U.S. Weather Bureau Res. Pap. 45, 58 pp.
Rae, W., Oxley, L., Reale, M., Brown, J., 2011. The emipirical properties of some
popular estimators of long memory processes, in 19th International Congress on
Modelling and Simulation, Perth, Australia, 2011.
52
Rallo, R., FerreGine, J., Arenas, A., Giralt, F., 2002. Neural virtual sensor for the
inferential prediction of product quality form process variables. Computers and
Chemical Engineering 26 (12), 173554.
Rampino, M.R., Self, S., 1982. Historic eruptions of Tambora (1815), Krakatau (1883),
and Agung (1963), their stratospheric aerosols, and climatic impact. Quaternary
Res. 18, 127143. doi: 10.1016/00335894(82)900655.
Randhir, T., Ekness, P., 2009. Urbanization on watershed habitat potential: A
multivariate assessment of thresholds and interactions. Ecohydrology 2, 88101.
Riebsame, W. E., Changnon, S.A., Karl, T.R., 1991. Drought and natural resources
management in the United States: Impacts and implications of the 198789 drought.
Westview Press New York, New York, 174 pp.
Ritter, H. Schulten, K., 1986. On the stationary state of Kohonen's selforganizing
sensory mapping. Biological Cybernetics 54 (1), 99106.
Said, S., Dickey, D.A., 1984. Testing for unit roots in autoregressivemoving average
models of unknown order. Biometrika 71, 599607.
Schwing, F.B., Murphree, T., Green, P.M., 2002. The Northern Oscillation Index (NOI):
A new climate index for the northeast Pacific. Progress in Oceanography 53, 115
139. http://www.esrl.noaa.gov/psd/gcos_wgsp/Timeseries/NAO/
Smith, T.M., Reynolds, R.W., 2003. Extended reconstruction of global sea surface
temperatures based on COADS data (18541997). Journal of Climate 16, 1495
1510.
Turcotte, D.L., 1997. Fractals and chaos in geology and geophysics, New York:
Cambridge University Press.
53
Ubilava, D. Helmers, C.G., 2013. Forecasting ENSO with a smooth transition
autoregressive model, Environmental Modelling & Software 40, 181190. ISSN
13648152, 10.1016/j.envsoft.2012.09.008.
Vesanto J., Alhoniemi E., 2000. Clustering of the SelfOrganized Map. IEEE
Transactions on Neural Networks 11(3), 586600.
Wang, S., 2003. Application of selforganizing maps for data mining with incomplete
data sets. Neural Computing and Applications 12, 4248.
Woodhouse, C.A., Overpeck, J.T., 1998. 2000 years of drought variability in the central
United States. Bulletin of the American Meteorological Society 79 (12), 26932714.
54

PAGE 1
FORECASTING CONDITIO NAL CLIMATE CHANGE USING A HYBRI D APPROACH By AKBAR AKBARI ESFAHAN I B.S. University of Colorado Denver 2010 A thesis submitted to the Faculty of the Graduate School of the University of Colorado in partial fulfillment of the requirements for the degree of Masters of Integrated Sciences Integrated Sciences 2013
PAGE 2
ii This thesis for the Masters of Integerated Science degree by Akbar Akbari Esfahani has been approved for the Integrated Sciences Program By Weldon Lodwick Chair Loren Cobb Michael J. Friedel September 23, 2013
PAGE 3
iii Akbari Esfahani, Akbar ( MSIS, Integrated Sciences ) M odeling Climate Change in the Southwestern United States Thesis directed by Professor Weldon Lodwick. ABSTRACT A novel approach is proposed to forecast the likelihood of climate change across spatial landscape gradients Th is hybrid approach involves reconstructing past precipi tation and temperature using the self organizing map technique; determining quantile trends in the climate change variables by quantile regression modeling; and computing conditional forecasts of climate change variables based on self similarity in quantile trends using the fractionally differenced auto regressive integrated moving average technique. The proposed modeling approach is applied to states (Arizona, California, Colorado, Nevada, New Mexico, and Utah) in the southwestern US, where conditional forecasts of climate change variables are evaluate d against recent (2012) observations, evaluate d at a future time period (2030), and evaluated as future trends (2009 2059) These results have broad economic, polit ical, and social implications because they quantify uncertainty in climate change forecasts affecting various sectors of society. Another benefit of t he proposed hybrid approach is that it can be extended to any spatiotemporal scale providing self similari ty exists. The form and content of this abstract are approved I recommend its publication. Approved: Weldon Lodwick
PAGE 4
iv DEDICATION I dedicate this work to the tireless support of my loving wife Naimeh and to my parents who made the biggest sacrifice by leaving behind their family so that we could have a better life.
PAGE 5
v ACKNOWLEDGMENTS I would like to thank Michael J. Friedel, my teacher, mentor and friend. For without his tireles s push, this research would have not been possible
PAGE 6
vi TABLE OF CONTENTS C hapter Table ................................ ................................ ................................ ................................ vii Figure ................................ ................................ ................................ ............................... viii 1. Introduction ................................ ................................ ................................ .................. 1 2. Methodology ................................ ................................ ................................ ................ 4 2.1 Reconstructing climate change variables ................................ ......................... 4 2.2 Quantile trends in climate change variables ................................ ..................... 7 2.3 Forecasts of climate change variables ................................ .............................. 8 2.3.1 Test of determinism ................................ ................................ ................. 8 2.3.2 Determining the fractal dimension ................................ ......................... 10 2.3.3 Long memory process ................................ ................................ ........... 11 2.3.4 Fractionally differenced auto regressive integrated moving average model ................................ ................................ ................................ ............... 12 2.4 Software for reconstruction, quantile trend s, and forecasting and empirical cumulative distribution functions ................................ ................................ ......... 13 3. R esults ................................ ................................ ................................ ........................ 14 3.1 Reconstructing climate change variables ................................ ....................... 14 3.2 Quantile trends in climate change variables ................................ ................... 20 3.3 Forecasts of climate change variables ................................ ............................ 21 3.3.1 A look at the past ................................ ................................ ................... 21 3.3.2 A look at the present ................................ ................................ .............. 23 3.3.2 A look at the future ................................ ................................ ................ 30 4. Conclusions ................................ ................................ ................................ ................ 46 References ................................ ................................ ................................ ......................... 47
PAGE 7
vii TABLES T able 3.1 Summary of long memory Hurst parameters for precipitation and temperature in Arizona, California, Colorado, New Mexico, Nevada, and Utah. ................................ .... 22 3.2 Forecast trends (2000 2050) for temperature, precipita tion, and drought in Arizona, California, Colorado, New Mexico, Nevada, and Utah. Results are for 0.25, 0.50, and 0.75 quantiles. ................................ ................................ ................................ ................... 43
PAGE 8
viii FIGURES Figure 2.1 Schematic depicting the hybrid modeling framework used to forecast climate change. ................................ ................................ ................................ ................................ 4 3.1 Quantile (upper panel) modeling of Palmer Drought Severity Index (PDSI) measurements (0 2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0.50 (gr een), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue) are indicated by double arrows. ................................ ................................ .............................. 16 3.2 Quantile (upper panel) modeling of reconstructed temperature measurements (0 2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue) are indicated by double arrows. ................................ ................................ ................................ ................................ ........... 17 3.3 Quantile (upper panel) modeling of reconstructed precipitation measurements (0 2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0 .50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue) are indicated by arrows. ....... 19 3.4 Comparison of conditional forecasts and measurement observations for 2012 Palmer Drought Severity index (PDSI): (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted (thin) and observed (open) PDSI value. ................................ .. 26 3.5 Comparison of conditional forecasts and measurement observations for 2012 temperature: (a) Cali fornia, (b) Colorado, and (c) Nevada. Arrows denote median forecasted (thin) and observed (open) temperature value. ................................ ................ 28 3.6 Comparison of conditional forecasts and measurement observations for 2012 precipitation: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted (thin) and observed (open) precipitation value. ................................ ............... 29 3.7 Conditional forecasts and histogram (lower panel) for Palmer Drought Severity Index in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 203 0 (thin) and observed 2012 (open) Palmer Drought Severity Index values. ................................ ................................ ................................ ............................... 32 3.8 Conditional forecasts and histogram (lower panel) for temperature in th e year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030 (thin) and observed 2012 (open) temperature values. ................................ ....................... 36 3.9 Conditional forecasts and histogram (lower panel) for precipitation in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030 (thin) and observed 2012 (open) precipitation values. ................................ ...................... 40
PAGE 9
1 1. Introduction People benefit from a multitude of resources and processes supplied by natural ecosystems (Randhir and Ekness, 2009). These benefits include water resources suitable for supporting various sectors of society such as agriculture, construction, daily living, energy, fishing, forestry, manufacturing, public health, recreation, transportation and overall economic development that maintains life systems prompting sustainability Climate change is frequently cited as one external driver of ecosystems (Furnis, 2010). Because climate is temporally and spatially dependent, change at a global scale differs from regional or local sc ales (Friedel, 2012 b ). One reason for spatial differences is the superposition of large scale climate patterns due to atmospheric and oceanic teleconnections (Schwing et al., 2002). Climate change also differs across temporal scales over which there are va riations in amplitude, gradient, and duration (Woodhouse and Overpeck, 1998; Mann et al., 2008; Mann et al., 2009; Friedel, 2011 a Friedel, 2012 a ). In many studies, the duration of climate change is considered short term (years to decades) variability (Woo dhouse and Overpeck, 1998, Ubilava and Helmers, 2013 ). Short term climate variability is often attributed to oscillations in the sea surface temperature (SST) that alter ocean currents and overlying air pressure resulting in a redistribution of temperature and precipitation (Smith and Reynolds, 2003). The El Nio Southern Oscillation (ENSO) is considered the strongest short term periodic fluctuation (2 7 years) with a rise (El Nio) or decrease (La Nia) of SST in the equatorial Pacific Ocean (Blade et al., 2008, Ubilava and Helmers, 2013). The influence of this teleconnection is not uniform in the United States, and E NSO events can affect things like water supply, water quality, riparian habitat, power generation, and range productivity. Related drought
PAGE 10
2 con sequences often include crop failure, debris flows, insect infestations, pestilence, violent conflict, wildfires, and disruptions to economic and social activities (Riebsame et al., 1991). Long term climate variability (hundreds to thousands of years) is o ften attributed to alterations in geologic and extraterrestrial processes, such as volcanic aerosols (Rampino and Self, 1982) and solar activity (Gray et al., 2010). Long term climate change reconstructions provide insight on past surface temperature and d rought variability at timescales crossing centuries or millennia (Mann et al., 2009; Friedel, 2012 a ). Should natural or anthropogenic forcing influence the frequency or intensity of climate change, there is an increased likelihood for future ENSO hazards p lacing national and global security at risk (Riebsame et al., 1991). For these reasons, climate change forecasts could benefit many sectors of society, but the scale dependent complexities render it a challenge using traditional process based models. Speci fically, climate forcing is known to interact with ecosystems characterized by coupled, nonlinear, and multivariate processes. Data associated with these ecosystems are typically sparsely populated ranging spatially from local (1000s km 2 ) to global and tem porally from immediate (1 10s years) to long term (100s to 1000 s years). One critical issue is the lack of essential calibration data that results in large inaccuracies (Loke et al., 1999). Other nonprocess based modeling efforts include linear time series models, such as the Autoregressive (AR) or Autoregressive I ntegrated M oving A verage (ARIMA) models ( Said and Dickey, 1984; Cowpertwait and Metcalfe, 2009) These traditional linear time series modeling schemes are too rigid with respect to detecting unexp ected features like the onset of nonlinear trends, or patterns restricted to sub samples of a data set.
PAGE 11
3 One alternative paradigm is to use a hybrid approach in which soft computing provides data for subsequent use in traditional numerical or empirical mod eling Some applications of the hybrid modeling approach are in rainfall runoff (Jain and Kumar, 2007), debris flow (Friedel, 2011 b ), mineral resource (Friedel, 2012 b ), and unexploded ordinance (Friedel et. al., 201 2 ). In this study, t he goal was to evaluate the efficacy of hybrid modeling for forecasting climate change over spatial landscape gradient s in the southwestern United States The objectives were to: (1) evaluate independent observations against the 2012 forecast of temperature, precipit ation, and drought across California, Colorado, and Nevada ; (2) evaluate future forecast trends over 50 years of probable temperature, precipitation, and drought across Arizona, California, Colorado, New Mexico, Nevada, and Utah ; and (3) to evaluate the pe rformance characteristics of models that are generated ( Bennett et al 2013) not just enough to evaluate the performance of a model, but it is a necessity to account for model uncertainty. It is to this end that t his study introduces a new modeling paradigm to account for uncertainty in time series forecasting. This study extends the work of Friedel (2011 a ) who sought to reconstruct 2,000 years of past temperature and precipitation for the south central and southwe stern United States ; Esfahani and Friedel ( 2010) who identified a long memory process in reconstructed climate variables ; Friedel (2012 a ) who used quantile regression to quantify uncertainty in global reconstructions of past temperature and precipitation and Caballero et al. (2002) and Nunes et al. (2011) who used a fractal approach to investigate the long memory process associated with temperature
PAGE 12
4 2. Methodology In modeling climate change variables, the hybrid approach relies on four computational steps: reconstruction, trends, forecasts, and uncertainty (Fig. 2. 1). Each of these steps is briefly described in the following sections. Fig 2 1 Schematic depicting the hybrid modeling framework used to forecas t climate change 2.1 Reconstructing climate change variables The reconstruction of climate change variables (temperature and precipitation) follows the approach described by Friedel (2012 a ). In that approach, a self organizing map (SOM) technique was used to project input data to a discrete lattice of competitive
PAGE 13
5 neurons (Kohonen, 2001), and estimates past climate change variables by minimizing topological error vectors (Fraser and Dickson, 2007). The process of projecting data is essentially a data c ompression technique (Hastie et al., 2002) for which the success of topology preservation was analyzed based on the quantization error E(G, X) given by : (1) where w i are weight vectors assigned to a fixed number of N neurons in the map grid G, x j are the M input data vectors (economic mineral resource variables), h i,I is a neighborhood function, is the Euclidian norm, and I is the best matching unit (BMU) vector. Implementation of the SOM learning method is based on the stochastic gradient described by Kohonen (2001). It consists of a two step process that is performed each time an input pattern is presented to the map: competition to determine the BMU and cooperat ive learning (spreading information contained in the current input vector across the map). At the beginning of the unsupervised training phase, the weight vectors are initialized to small random numbers. The input data vectors are presented to the map grid in a random fashion to generate data clusters without introducing bias for a specific class. In the first step, the BMU with map coordinates ( I i I j ) is determined as the grid neuron, whose weight vector is the closest to the input given by : (2) where ArgMin is the minimum distance defining the centr al position of the neighborhood function. The neighborhood function h i, I is chosen t o be Gaussian function given by :
PAGE 14
6 (3) w here r i r 1  corresponds to the distance between map neuron r i and BMU in the map grid, and (n) defines the width of the neighborhood function, a monotonically decreasing function of the iteration (also called epoch) number n. In the second step, a weight update is determined which is a function of the distance to the current BMU, as expressed through the neighborhood function h i,I (n). The weights are gradually adjusted according to: (4) where is a scalar value called the learning rate bounded on the interval [0, 1]. The B MU ensures that the largest weight correction ( h i,I (n) =1) is adjusted in the direction of the input vector. The association effect takes place at the neighboring nodes but to a l esser degree because of the Gaussian shape. This adaption procedure stretches the w eight vectors of the BMU and its topological neighbors towards the input vector. Presenting similar input vectors to the map provides further activations in the same neig hborhood and thereby tends to produce clustering of data in the feature space. Association between neurons decreases during the learning process (the width of the neighborhood function (n) is forced to decrease with n preserving large clusters of data whi le enabling the separation of clusters that are closely spaced). Ultimately, this training process results in a topology where similarities among data patterns are mapped into similar weights of the neighboring neurons, and the asymptotic local density of the weights approach that of the training set (Ritter and Schulten, 1986).
PAGE 15
7 Cross validation (Efron and Tibshirani, 1993) is conducted to ensure the SOM provides unbiased estimates of climate change variables. In this case, known data values are estimated based on distances among the available model vectors (Wang, 2003; Kalteh and Berndtsson 2008). In the traditional approach, estimates of values are taken directly from the prototype vectors of the best matching units (Fessant and Midenet, 2002; Wang, 2003) Often times certain training data sets result in biased estimates (Dickson and Giblin, 2007; Malek, et al., 2008) requiring a modified scheme that incorporates bootstrapping (Breiman, 1996), ensemble average (Rallo et al., 2002), or nearest neighbor (Mal ek et al., 2008). This study uses an alternative iterative estimation scheme that minimizes the topological error vector (Fessant and Midenet, 2002). The estimation of past climate change values for all variables is done simultaneously and referred to here as the reconstruction. For more details about SOM training and estimation, the reader also is referred to (Kohonen, 2001; Vesanto, and Alhoniemi, 2000). 2.2 Quantile trends in climate change variables The determination of quantile trends in the reconstruc ted climate change variables follows the approach described by Friedel ( 2012) One advantage of this approach is its flexibility in modeling data with conditional functions that may have systematic differences in dispersion, tail behavior, and other covari ate features (Koenker, 2005). In adopting this approach the quantile regression approach is adopted thereby reduc ing the d dimensional nonparameteric regression problems to a series of additive univariate problems given by :
PAGE 16
8 (5) where min is the minimization operator and estimated coefficients at distinct is an n dimensional vector with elements in the j, j+1, and j+2 positions and zeros elsewhere ; is the quantile regression loss function; is a factor that controls the degree of smoothing; T is the transpose; are response observation values associated with random variable X; and are response observation values associated with random var iable Y. Individual quantile curves can then be specified as a linear b spline of the form : (6) where is the th conditional quantile function, are the regression quantiles; are regression coefficients; denote the basis function of the spline; and is the quantile (0 < < 1). Once the knot positions of the spline have been selected, such models are then linear in parameters and therefore can be estimated (Koen ker, 2005). 2.3 Forecasts of climate change variables 2.3. 1 Test of determinism Before fractal techniques can be used to analyze the reconstructed climate change data, the determinism of a fractal time series must be established T he data must be identified as deterministic or statistical (Kaplan, 1994; Turcotte, 1997). The other assumption to be satisfied is that stationarity exists for long memory process data (Beran,
PAGE 17
9 1994). The determinism test developed by Kaplan (1994), allows the determination as to whether a time series has a deterministic or stochastic structure. Turcotte (1997) defines the set as a deterministic fractal set, if the set is scale invariant at all scales; and a statistical fractal set, if the set is different a t different scales but the differences do not allow the scale to be determined. This is an important distinction because it determines the correct way to calculate the fractal dimensions (as defined by equation (10)) for a time series. Mathematically, the delta epsilon test of continuity (Kaplan, 1994) is applied to determine the deterministic structure of a time series. To do so, the continuity test is applied to orbits comprising the phase space topology created by time delayed embedding of the original set. This process is facilitated by generating an ensemble set of surrogate time series using a bootstrapping approach ( Constantine et al., 2010 ) Next, the phase space statistic, called the E statistic (Kaplan, 1994), is calculated for the time delayed embedding of the original time series as well as the ensemble of surrogate data. The structure is then judged by comparing the E statistics of the original set to the E statistic of the surrogate set. A separation between the two statistics implies the exi stence of a deterministic structure, and the converse implies the set is a realization of a random process and thus has a stochastic structure. According to Kaplan (1994) the E Statistic is defined as: (7) (8) (9)
PAGE 18
10 where is the Euclidean distance between phase space points z(j) and z(k), and is the corresponding separation distance between the points at a time points in the future along their respective orbits. The variable is the orbital lag and the future points are the images of the original pair. The increment is the width of a specified Euclidean bin size. Given the distance is used to ident ify the proper bin in which to store the image distance and the average of each bin forms the statistic. Finally, the E statistic is formed by calculating the cumulative sum of the statistic ( Constantine and Percival, 2011) 2. 3.2 Determining th e fractal dim ension Fractal sets are characterized by their dimensions ( Mandelbrot, 1967; Mandelbrot and van Ness, 1968) Fractal dimensions are mathematically different from topological dimensions in that while topological dimensions are strictly constructed from integers, fractal dimensions can be fractional. While there are an infinit e number of fractal dimensions, t he fractal information dimension is important because it describes the entropy of a data set. In general, the fractal dimension is d escribed by: (10) where q refers to the dimension of interest. By embedding the dataset in an n dimensional grid with cells having sides of size r, we compute the frequency for which a data point falls into the i th cell, p of I. The information dimension is described by As the information dimension approaches 1, the numerator of equation (4) changes to the Barbara,
PAGE 19
11 1999) Un der these conditions, the fractal dimension of a stochastic self similar data set simplifies to: (11) Without the partial differential equation requirement the problem becomes computationally tractable and can be calculated easily, w hile there are several ways to calculate the fractal dimensions, f or the purposes of this paper, the box counting technique is used ( De Pison et al., 2008) 2. 3.3 Long memory process The existence of a long memory process was first explored by Hurst (1951) when trying to find a solution on how to regularize flow of the Nile River. He observed that long periods of high flow levels were followed by long periods of low flow levels indicating the existence of a long memory process. In general, l ong m emory p rocess describes the process where by a time series shows that past values have a strong effect on present or future values. Mandelbrot and van Ness (1968) introduced the Hurst parameter (H) to describe the long term memory of a time series process. The Hu rst parameter is given by : (12) where is the long memory parameter with the range of thus, H has a range of The closer H is to 1, the more persistent the time series is considered and at values less than or equal to a long me mory process does not exist. Currently there are several methods to estimate the Hurst parameter; however, Rae et al. (2011) showed
PAGE 20
12 1994). 2. 3.4 Fractionally difference d auto regressive integrated moving average model If climate change variables have a long memory process then the fractionally differenced a uto regressive integrated moving average model (FARIMA) can be applied and forecasts generated. The FARIMA model is similar to the Box Jenkins ARIMA model ( Cowpertwait and Metcalfe, 2009) except that the integrated part of the ARIMA model can be a fractional number defined (Beran, 1994) as: (13) This expression underscores the need for a good estimation of the Hurst parameter when dealing with climate time series data. If the climate change data exhibit s long memory process and stationarity, the FARIMA model assumptions are satisfied and can be us ed to generate forecasts following three steps. First, the model parameters for each time series are determined using the maximum likelihood estimators ( Fraley et al., 2011) Second, the model is fitted and forecasts made using the estimated parameters ( Hy ndman, 2011) Third, the forecasts are validated using various statistical measures (Hyndman and Koehler, 2006). 2.3.5 Empirical cumulative distribution function The use of quantile regression facilitates quantization of the prediction uncertainty by const ructing empirical cumulative distribution functions (ECDF) at any year of a forecast. The empirical cumulative distribution function (ECDF) is determined by modeling the collection of quantile forecast results. The ECDF is a step function F n :
PAGE 21
13 (14) where i is the number of tied observations at that value (missing values are ignored). For observations and is the fraction of observations less or equal to t. 2.4 Software for reconstruction, quantile trends and forecasting and empirical cumulative distribution functions The data mining, reconstruction, and analysis are carried out using the SiroSOM (CSIRO Exploration & Mining, 2008) graphical user interface (GUI). This GUI provides an interface between data sets and functions in the freely available SOM Toolbox (Adaptive Informatics Research Center, 2010). Quantile modeling is conducted using the quantreg, splines, and stats packages (Hornik, 2011); fractal transformation and fractal dimension calculations a re conducted using the fdim package (Martinez de Pison, Ascacibar et al. 2008). Fractionally differenced auto regressive integrated moving average modeling is conducted using the forecast and fractal package (Constantine and Percival, 2011, Hyndman 2011) a nd the empirical cumulative distribution function modeling is conducted using the stats package (Hornik, 2011); freely available at http://www.r project.org/ in the R toolbox
PAGE 22
14 3. R esults 3.1 Reconstructing climate change variables Annual temperature and precipitation values were reconstructed across a gradient of modern climate zones: Arizona (Desert), California (Mediterranean), Colorado (Semiarid to Alpine), Nevada (Semiarid to Arid), New Mexico (Semiarid), and Utah (Semiarid). Th e simultaneous reconstruction of annual temperature and precipitation was done based on the self organized nonlinear data vector relations among approximately 2,000 years (0 to 2009 AD) of reconstructed warm season (average of June, July, and August) Palme r Drought Severity Index (PDSI) data (Cook et al., 2004), 114 years (1895 2009) of annual state precipitation (accumulation over January through December) and temperature (average of January through December) data (National Climatic Data Center, 2010), an d other related tropical and extratropical measurements described by Friedel (2012a). The temperature and precipitation are standard climate variables, whereas the PDSI defines annual dry, neutral, and wet periods based on tree ring information (Palmer, 19 65; Cook et al., 2004). For a comprehensive review of the reconstruction and validation of these data, the reader is referred to the cited references. In this study, the reconstructions were verified against independent precipitation and temperature data f or the years: 1896, 1900, 1911, 1919, 1923, 1935, 1940, 1952, 1960, 1966, 1968 (La Nia), 1986, 1993, 1998 (La Nia), and 2005 (El Nio) using split and cross validation (leave one out) approaches (Bennett et al, 2013). The Spearman Rho correlation among observed and reconstructed values was greater than 95% with a p value of 0.001. The magnitude of climate variability over the past 2,000 years is visually apparent when inspecting plots for PDSI (Fig. 3.1 a, b, and c), temperature (Fig. 3 .2 a, b,
PAGE 23
15 a nd c), a nd precipitation (Fig. 3.3 a, b, and c), where (a) is California, (b) is Colorado, and (c) is Nevada. (a) California (b) Colorado
PAGE 24
16 (c) Nevada Fig. 3 1 Quantile (upper panel) and histogram (lower panel) modeling of Palmer Drought Severity Index (PDSI) measurements (0 2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the Medieval Warming Pe riod (red) and Little Ice Age (blue) are indicated by double arrows. (a) California
PAGE 25
17 (b) Colorado (c) Nevada Fig. 3 2 Quantile (upper panel) and histogram (lower panel) modeling of reconstructed temperature measurements (0 2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue) are indicated by double arrows.
PAGE 26
18 (a) California (b) Colorado
PAGE 27
19 (c) Nevada Fig. 3 3 Quantile (upper panel) and histogram (lower panel) modeling of reconstructed precipitation measurements (0 2009): (a) California, (b) Colorado, and (c) Nevada. The trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles. The approximate timing of the Medieval Warming Period (red) and Little Ice Age (blue) are indicated by arrows. In general, the plots reveal that all states in the southwest experienced past conditions ranging from extreme drought to extremely moist. It is interesting to note that occ urrence in California (62.5%), Colorado (58.5%), and Nevada (55%). This suggests that while there is the likelihood for additional extreme events, there is a greater likelihood for current drought conditions to be eventually mitigated. Regarding temperatur e, the respective mode values in California, Colorado, and Nevada are about 14.8 (range from 14 to 16 ), 6.75 (range from 5.5 to 9.5 ), and 9.75 (range from 8 to 12 ) W hereas the respective mode values for precipitation in California, Color ado, and Nevada are about 425 mm (range from 150 mm to 1050 mm),
PAGE 28
20 350 mm (range from 230 mm to 540 mm), and 190 mm (range from 120 mm to 350 mm). Both the temperature and precipitation for these states reflect their association with the modern climate gradi ent. 3.2 Quantile trends in climate change variables Q uantile regression modeling is applied to the past (0 2009) PDSI (Fig. 3.1 a, b, &c), temperature (Fig. 3 .2 a, b, &c) and precipitation (Fig. 3.3 a, b, &c) data. The annual quantiles trends are presented for climate change data in a) California, b) Colorado, and c) Nevada. The decadal trends reflect 0.05 (red), 0.25 (yellow), 0.50 (green), 0.75 (blue), 0.95 (purple) quantiles determined using 200 b spline degrees of freedom. In general, t he quantile trends reveal that the long term regional climate was interrupted by short term changes. T he so called Medieval Warm Period (~900 to ~1250) and Little Ice Age (~1 3 00 to ~1850) appear as primary lower frequency disruptio ns over the last two millennia in California and Nevada, and secondary (muted) disruptions in Colorado (Crowley and Lowery, 2000; Mann, 2002). The approximate timing of the Medieval Warming (horizontal red arrow) coincides with a decrease in PDSI (drier co nditions), increase in temperature, and decrease in precipitation; whereas the Little Ice Age Periods (horizontal blue arrow) coincides with an increase in PDSI (wetter conditions), decrease in temperature, and increase in precipitation. These findings are attributed to strong ENSO teleconnection s with California and Nevada, but mixed ENSO signals in Colorado ( it is a region between El Nio and La Nia latitudes ). In computing forecasts for future climate in southwestern states, additional quantile trends ( 0.05, 0.10, 0.20, 0.30, 0.40, 0.50, 0.60, 0.70, 0.80, 0.90, and 0.95 ) are used to increase resolution when modeling the empirical cumulative distribution functions.
PAGE 29
21 3.3 Forecasts of climate change variables 3.3 1 A look at the past To evaluate past climate change data for long memory process ( Barbara, 1999) the PDSI, temperature, and precipitation are transformed to an equivalent fractal information dimension The transformation is applied to 20 year intervals resulting in time series plots comprising 100 lags. The various methods used to calculate the long memory Hurst parameter (H) are summariz ed in table 3.1 Based on these H calculations, the precipitation appears to have the strongest long memory process, whereas temperature and PDSI have a comparative ly weak er long memory process. This finding implies that the autocorrelation of precipitation is stronger than temperature ; that is, the precipitation persists over a longer number of years than temperature Another finding is that the arithmetic calculation of H falls within the bounds of the Whittle approximation. One exception may be calculations based on the spectral regression method. For example, application of this method to precipitation results in an H parameter value that is outside the bounds of the Whittle approximation. This suggests that the spectral regression method might not provide reasonable approximations for comparatively strong long memory processes. It is also interesting to note that the two driest states, Nevada and Arizona, have the lowest H values for precipitation. This suggests that wet cycles in these states are of shorter duration than other southwestern states. By contrast, Colorado has the second shortest temperature cycle and Ca lifornia the longest temperature cycle. At the same time, Colorado has the longest drought cycle and California has the shortest drought cycle. These findings suggest that Colorado is more likely to experience extended drought cycles at different temperatu res, whereas
PAGE 30
22 California is more likely to experience longer cycles of high temperatures that do not result in droughts which can be explained by the fact that California borders the Pacific Ocean Table 3 1 Summary of l ong memory Hurst parameter s for precipitation and temperature in Arizona, California, Colorado, New Mexico, Nevada, and Utah. Precipitation H a calculation H interval Whittle Method H via Spectral Regression 2.50% 97.50% AZ 0.948 0.850 1.171 0.69 3 CA 0.97 8 0.86 9 1.234 0.701 CO 0.97 9 0.871 1.190 0.74 2 NM 0.98 3 0.9 30 1.368 0.742 NV 0.96 4 0.83 2 1.152 0.69 3 UT 0.982 0.934 1.31 2 0.999 Temp erature H a calculation H interval Whittle Method H via Spectral Regression 2.50% 97.50% AZ 0.60 6 0.569 0.651 0.616 CA 0.662 0.54 6 0.7 90 0.59 7 CO 0.57 9 0.53 3 0.62 8 0.59 2 NM 0.60 2 0.521 0.73 8 0. 600 NV 0.560 0.52 1 0.61 1 0.570 UT 0.60 5 0.64 7 0.895 0.60 1 Palmer Drought Severity Index H a calculation H interval Whittle Method H via Spectral Regression 2.50% 97.50% AZ 0.63 9 0.5 90 0.696 0.64 7 CA 0.554 0.507 0.60 4 0.53 8 CO 0.651 0.596 0.707 0.68 2 NM 0.633 0.57 4 0.692 0.649 NV 0.630 0.572 0.68 8 0.64 6 UT 0.59 8 0.533 0.66 4 0.638 3.3.2 A look at the present The previous section established that the southwest climate change reconstructions are stationary and stochastically self similar. Given th ese facts, the
PAGE 31
23 FARIMA model is applied to individual quantile trends for forecasting of climate change variables. In all cases, t he fitted FARIMA models are characterized as one of the following types: from AR(0) MA(0) to AR(0) MA(6). Using the fitted FARIMA models, forecasts are generated over the range of quantiles from which ECDFs are then computed for each variable. T he FARIMA model performance is evaluated by comparing the median 2012 state forecasts of PDSI (Fig. 3.4 ), temperature (Fig. 3.5 ), and precipitation (Fig. 3.6 ) (California, Colorado, and Nevada respectively) to state observations ( National O ceanic and Atmospheric Administration, 2013 ). In this comparison, the PDSI reflects values averaged over the period of June July August, the temperature represents values averaged over the period of January December, and the precipitation reflect values ac cumulated precipitation over the period of January December. When comparing the 2012 climate change observations to forecasts, it is important to note the differences in their spatiotemporal representation. Specifically, the observed PDSI and reconstructed climate variables are values associated with a grid location (influenced by a 250 km x 250 km region) defined by Cook et al. (2004), whereas the observed climate variables represent values averaged over climate divisions with larger areas defined by the N ational Oceanic and Atmospheric Administration ( 2013 ) The respective 2012 state observations and forecasts of PDSI, temperature, and precipitation are presented in Figs. 3.4 3.6 Reviewing the plots for C alifornia Colorado and Nevada reveals that all ob servations plot within the probable forecast limits supporting the usefulness of the hybrid approach. In addition the range of probable climate variable differs among states reflecting the modern climate gradient and
PAGE 32
24 underscoring the stationarity of conditions over the past 2000 years. For example, the respectiv e forecast range of PDSI (Fig. 3.4 ) for California, Colorado, and Nevada is from 3.50 to 0.5, from 4.25 to 2.25, and from 5.0 to 1.0. The respective forecast range of temperature (Fig. 3.5 ) for California, Colorado, and Nevada is from 15.1 to 15.8 from 6.65 to 8.0 and from 9.6 to 10.5 The respective forecas t range of precipitation (Fig. 3.6 ) for California, Colorado, and Nevada is from 200 mm to 625 mm, from 340 mm to 440 mm, a nd from 170 mm to 260 mm. Comparing the median forecast (indicated by vertical arrows associated with the 0.5 quantile) to observed climate variables averaged over different climate divisions reveals reasonable correspondence. For example, the respective median forecast and average d observations of PDSI in California (Fig. 3.4 a) are about 2.50 and 2 .20 (climate divisions 4 and 6). T he respective median forecast and averaged observation of PDSI in Colorado (Fig. 3.4 b) are about 0.75 and 2.71 (climate divisions 4 and 1); and the r espective median forecast and averaged observation of PDSI in Nevada (Fig. 3.4 c) are about 2.15 and 2.17 (climate divisions 3 and 4) Similarly, the respective median forecast and averaged observation of temperature in California (Fig. 3.5 a) are about 15 .44 and 15.2 (climate divisions 4 and 6) T he respective median forecast and averaged observation of temperature in Colorado (Fig. 3.5 b) are about 7.42 and 7.78 (climate divisions 4 and 1); and the respective median forecast and averaged observatio n of temperature in Nevada (Fig. 3. 5 c) are about 10.05 and 10.6 (climate divisions 3 and 4). The respective median forecast and averaged observation of pre cipitation in California (Fig. 3.6 a) are about 450 mm and 440 mm (climate divisions 4 and 6). T he respective median forecast and averaged observation of precipitation in Colorado (Fig. 3.6 b) are
PAGE 33
25 about 403 mm and 392 mm (climate divisions 4 and 1); and the respective median forecast and averaged observation of precipitation in Nevada (Fig. 3.6 c) are ab out 202 mm and 188.5 mm (climate divisions 3 and 4). Comparison of median forecasts to values for individual divisions and all state divisions reveals the heterogeneous nature of climate variables within each state; that is, the observations across climate regions in California are mostly Mediterranean; in Colorado, they range from s emiarid to a lpine ; in Nevada they range from s emiarid to a rid This finding suggests that future improvements, potentially spanning five to 15 years, to forecasting may be achieved by introducing PDSI data associated with grid nodes across all state climate divisions. (a) California
PAGE 34
26 (b) Colorado (c) Nevada Fig. 3 4 Comparison of conditional forecasts and measurement observat ions for 2012 Palmer Drought Severity index (PDSI): (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted (thin) and observed (open) PDSI value.
PAGE 35
27 (a) California (b) Colorado
PAGE 36
28 (c) Nevada Fig. 3 5 Comparison of conditional forecasts and measurement observations for 2012 temperature: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted (thin) and observed (open) temperature value. (a) Californi a
PAGE 37
29 (b) Colorado (c) Nevada Fig. 3 6 Comparison of conditional forecasts and measurement observations for 2012 precipitation: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted (thin ) and observed (open) precipitation value.
PAGE 38
30 3.3 .2 A look at the future The conditional state 2030 forecasts are presented together with 2012 observations (average of climate divisions) for PDSI, temperature, and precipitation (Figs. 3.7 3. 9 ). For example, the conditional forecasts of PDSI in California, Colorado, and Nevada are presented in Fig. 3.7 (a) California
PAGE 39
31 (b) Colorado
PAGE 40
32 (c) Nevada Fig. 3 7 Conditional forecasts and histogram (lower panel) for Palmer Drought Severity Index in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030 (thin) and observed 2012 (open) Palmer Drought Severity Index values. In California (Fig. 3.7 a ), the median PDSI forecast value of about 1.42 is larger than the observed 2012 value of about 2.20 (average of climate divisions 4 and 6). This
PAGE 41
33 suggests the likelihood for a shift toward moister conditions; that is, the current moderate drought conditi on is likely to shift to a mid range condition within the next 8 years. In Colorado (Fig. 3.7 b), the median PDSI forecast value of about 0.40 is larger than the observed 2012 value (average of climate divisions 4 and 1) of about 3.50, also suggesting a s hift from severe drought to a mid range condition in the next 8 years. In Nevada (Fig. 3.7 c), the median PDSI forecast value of about 1.60 is larger than the observed 2012 value of about 2.17 (average of climate divisions 4 and 6), also suggesting the li kelihood for a shift toward moister conditions in the next 8 years. In addition to trends based on median values, there also is some probability that California could experience conditions ranging from severe drought to moderately moist in the next 8 years Colorado could experience conditions ranging from extreme drought to moderately moist in the next 8 years and Nevada could experience conditions ranging from severe drought to mid range in the next 8 years. Relative to Colorado, California tends toward t he moist and Nevada shifts toward drought conditions.
PAGE 42
34 (a) California
PAGE 43
35 (b) Colorado
PAGE 44
36 (c) Nevada Fig. 3 8 Conditional forecasts and histogram (lower panel) for temperature in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030 (thin) and observed 2012 (open) temperature values. The conditional state 2030 forecasts of temperature in California, Colorado, and Nevada are presented in Fig. 3.8 In California (Fig. 3.8 a), the median temperature
PAGE 45
37 forecast value of about 15.28 is about the same as the observed 2012 temperature (average of climate divisions 4 and 6) of about 15.2 This suggests the likelihood for a slight shift toward warmer conditions in the next few years. In Colorado (Fig. 3.8 b), the median temperature forecast value of about 7.3 is less than the 2012 temperature (average of climate divisions 4) of about 7.78 suggesting a shift toward cooler conditions in the next few years. In Nevada (Fig. 3.8 c), the median temperature forecast value of about 9.9 is smaller than the 2012 temperature (average of climate divisions 4 and 6) of about 10.6 suggesting the likelihood for cooler conditions in the next few years. In addition to trends based on medi an values, there also is some probability that California could experience temperatures ranging from about 14.7 to 15.8 Colorado could experience temperatures ranging from about 6.5 to 8.25 and Nevada could experience temperatures ranging from ab out 9.2 to 10.7 Within the next two to eight years, relative to Colorado, California and Nevada tend toward warmer conditions.
PAGE 46
38 (a) California
PAGE 47
39 (b) Colorado
PAGE 48
40 (a) Nevada Fig. 3 9 Conditional forecasts and histogram (lower panel) for precipitation in the year 2030: (a) California, (b) Colorado, and (c) Nevada. Arrows denote median forecasted 2030 (thin) and observed 2012 (open) precipitation values. The conditional state 2030 forecasts of precipitation in California, Colorado, and Nev ada are presented in Fig. 3.9 In California (Fig. 3.9 a), the median precipitation
PAGE 49
41 forecast value of about 490 mm is larger than the 2012 precipitation value (average of climate divisions 4 and 6) of about 440 mm. This suggests the likelihood for a future shift toward drier conditions in the next few years. In Colorado (Fig. 3.9 b), the median precipitation forecast value of about 390 mm is about the same as the 2012 precipitation value (average of climate divisions 4) of about 392 mm, suggesting no change i n precipitatio n conditions. In Nevada (Fig. 3.9 c), the median temperature forecast value of about 205 mm is larger than the historical median value of about 188.5 mm, suggesting the likelihood for future wetter conditions in roughly eight years. In additio n to trends based on median values, there also is some probability that California could experience precipitation amounts ranging from 275 mm to 725 mm, Colorado could experience precipitation amounts ranging from 300 mm to 470 mm, and Nevada could experie nce precipitation amounts ranging from 150 mm to 300 mm within the next eight years. Relative to Colorado, California tends toward wetter conditions and Nevada toward drier conditions. In addition to evaluating forecasts for the year 2030, an analysis of t rends is conducted over forecasts spanning a 50 year period ( from 20 10 to 20 60 ). A summary of findings is presented for the southwestern states in Table 2 These results indicate that the southwestern US is likely to experience decreasing temperatures that will likely continue over the next 25 years resulting in wetter and cooler conditions for the region The decay represents a decrease of 0.20 degrees Celsius over 20 years in temperature. The accuracy of the forecast results for temperatures are verified by comparing the root mean squared error (RMSE) ( Bennett et al 2013 ; Hyndman and Koehler, 2006 ) of the forecast of each
PAGE 50
42 state to the mean of the forecast of that state. For temperature, the range is from 2.3% to 8.2% of RMSE to the mean indicating a low variability. The PDSI index will increase in the same period indicating wetter conditions across the southwestern states. However, this increase is only by about 0.8 on the PDSI index, which suggests only a very moderate decrease of drought. The exception to this increase is California where the forecast suggests a constant level of PDSI suggesting that California has reached its equilibrium, the end of a long memory process. Since PDSI values are mainly negative, RMSE is not a good indicator of performanc e, since it relies on the absolute values, instead the mean error (ME) of the residuals is chosen here. A comparison between the ME of each state and the mean of the forecast of each state indicates variability of about 2% to 4%.
PAGE 51
43 Table 3 2 Forecast trends (2000 2050) for temperature, precipitation, and drought in Arizona, California, Colorado, New Mexico, Nevada, and Utah. Results are for 0.25, 0.50, and 0.75 quantiles. Palmer Drought Severity Index (PDSI) Forecast Accuracy Mean of Series Forecast Direction State ME 1 AZ 0.020 0. 54 increase CA 0.016 0. 24 constant CO 0.011 0. 51 increase NM 0.014 0. 61 increase NV 0.018 0. 70 increase UT 0.009 0.45 increase Increase = increase of 0.8 to 1.0 on PDSI scale over 10 years Decrease = decrease of 0.8 to 1.0 on PDSI scale over 10 years Since PDSI is measured from positive to negative, increase indicates wetter conditions Precipitation Forecast Accuracy Mean of Series Forecast Direction State RMSE 2 MAE 3 AZ 59.5 44.2 306.1 sharp increase CA 98.3 77.0 528.2 increase CO 45.7 35.7 387.7 decrease NM 55.6 41.9 333.8 increase NV 31.3 24.9 215.0 decrease UT 39.4 31.1 283.9 increase Increase = increase of 5 mm over 10 years Sharp Increase = increase of 20 mm over 10 years Decrease = decrease of 5 mm over 10 years Sharp Decrease = decrease of 20 mm over 10 years Temp erature Median of Series Forecast Accuracy Mean of Series Forecast Direction RMSE MAE AZ 0.6 2 0.51 15. 6 decrease CA 0.39 0.32 15. 1 decrease CO 0.59 0.47 7. 2 decrease NM 0.46 0.3 7 11.8 decrease NV 0.58 0.4 5 9. 7 decrease UT 0.60 0.50 9.0 decrease Decrease = decrease of 0.20 degrees over 20 years Increase = increase of 0.20 degrees over 20 years 1 ME mean error, since PDSI is mostly a negative value, MAE and RMSE do not apply here. 2 RMSE root mean squared error. 3 MAE mean absolute error.
PAGE 52
44 The forecast results for precipitation is slightly different than those of temperature and PDSI. The forecast accuracy is measured using the RMSE for each state, which indicates a variability of 11% to 19%. The higher variability can be explained by the fa ct that precipitation has much stronger long memory cycles (150 to 300 years) and a fifty year forecast is only part of the cycle. It is interesting to note that Nevada and Colorado will experience an increase in precipitation in the near future (next 5 to 10 years), however; the long term forecast horizon (10 to 20 years) calls for less precipitation. This indicates that the short forecast horizon for these 2 states is at the end of one cycle and the beginning of the next cycle.
PAGE 53
45 4. Conclusions The proposed hybrid modeling approach is useful for forecasting drought, temperature, and precipitation at the level of state climate divisions for a span of 1 to 15 years. This process requires the reconstruction of past climate variables that are stationar y and stoc hastically self similar with quantile regression modeling used to facilitate quantization of forecast uncertainty. The application to southwestern states provided reconstructions of past (0 2009) climate records exhibiting temporal paleoclimatic features such as the Medieval Warm Period and Little Ice Age. Independent p erformance testing using modern (2012) observations averaged over appropriate climate divisions demonstrated good correspondence to median forecasts. This finding supports the 50 years of f orecasting (2010 2060) to assist managers in formulating decisions, potential mitigating strategies, and policy associated with future, uncertain climate change Differences among median forecasts and observations from other climate divisions demonstrate t he heterogeneous nature of climate variables within each state. This finding supports future improvements in forecasting by introducing additional paleoclimatic data associated with grid nodes crossing all of the state climate divisions. Because t he propos ed hybrid approach can be extended to any spatiotemporal scale providing self similarity exists forecasting has the possibility to address economic, political, and social aspects affecting various sectors of society As with any data driven approach, the introduction of additional related information can be expected to further reduce the uncertainty. While the approach presented was geared toward climate change, the work here can be extended to any type of time series analysis. It re presents a superior method to
PAGE 54
46 traditional accuracy measures of forecasting. That is, traditional measures are tuned to account for errors in the model, however; the method presented here accounts for the uncertainty of the entire modelling process. The met hod can be applied to time series and spatial analysis in fields such as economics, ecology, and biology. Beyond the novel method presented for forecasting uncertainty, the fractal information dimension transformation of data can be used to find changes to time series. Current research in the area concetrates on the correlation of volcanic acitivities to climate shocks found in the fractal transformation of 2000 years of data presented in this research (Friedel and Akbari Esfahani, 2013).
PAGE 55
47 REFERENCES Akbari Esfahani, A., Friedel, M. J. 2010 The fractal nature of climate change 2000 years in retrospect H21G 1126 presented at 2010 Fall Meeting, AGU, San Francisco, Calif., 13 17 Dec. Akbari E sfahani, A. Friedel, M.J., 2012, A Fractional Forecast to the Climate of Southwestern United States, Mathematics of Planet Earth presented at 2012 Meeting, IUGG Conference on Mathematical Geophysics, Edinburgh, United Kingdom, 18 22 June. Adaptive Informa tics Research Center (AIRC) 2010. SOM Toolbox, Helsinki University of Technology, Laboratory of Computer and Information Science, Adaptive Informatics Research Center. http://www.cis.hut.fi/projects/somtoolbox/ Barbara, D., 1999. Chaotic mining: Knowledge discovery using the fractal dimension, George Mason University Information and Software Engineering Department, Fairfax, VA. Bennet, N.D., Croke, F.W., Guariso, G., Guillaume, J.H.A., Hamilton, S.H., Jakeman, A.J., Marsili Libelli, S., Newham, L.T.H., N orton, J.P., Perrin, ., Pierce, S.A., Robson, B., Seffelt, R., Voinov, A.A., Fath, B.D., Andressian, V., 2013, Characterising performance of environmental models, Environmental Modelling & Software 40 (2013), 1 20. Beran, J., 1994. Statistics for long mem ory process, Boca Ranton: CRC Press LLC. Blade, I., Newman, M., Alexander, M.A., Scott, J.D., 2008. The late fall extratropical response to ENSO: Sensitivity to coupling and convection in the tropical west Pacific. J. Climate 21(23), 6101 6118.
PAGE 56
48 Breiman, L. 1996. Bagging predictors, Machine Learning, 24, 123 140. Caballero, R., Jewson, S. Brix, A., 2002. Long memory in surface air temperature: Detection, modeling, and application to weather derivative valuation. Climate Research 21, 127 140. Constantine W., Percival D., 2010. fractal: Insightful Fractal Time Series Modeling and Analysis. R package version 1.0 3. http://CRAN.R project.org/package=fractal Cook, E. R., Woodhouse, A. Meko, D. M. Stahle, D. W., 2004. Long term aridity changes in the western United States, Science 306 (5698), 1015 1018. Constantine, W., Percival, D., 2011. Fractal: fractal time series modeling and analysis. R package version 1.1 1, http://CRAN.R project.org/package=fractal, 2011. Cowpertwait, P.S., Metcalfe, A.V., 2009. Intr oductory time series with R, New York: Springer. Crowley, T.J., Lowery, T.S., 2000. How Warm Was the Medieval Warm Period? AMBIO: A Journal of the Human Environment 29 (1), 51 54. CSIRO Exploration & Mining, 2008, SIROSOM software. http://www.csiro.au/science/SOM data analysis.html Dickson, B.L., Giblin, A.M., 2007. An evaluation of methods for imputation of missing trace element data in groundwaters. Geochemistry: Exploration, Environment, and Analysis 7 (2), 173 178. Efron B., Tibshirani R.J., 1993. An introduction to the bootstrap, in Monographs on Statistics and Applied Probability, 57, 436 pp, Chapman & Hall, London. Fessant, F., Midenet, S., 2002. Self organizing map for data imputation and correction in surveys. Neural Comp uting Applications 10, 300 310.
PAGE 57
49 Fraley, Leisch, F., Maechler, M., Reisen, V., Lemonte, A., 2011. Fracdiff: Fractionally differenced ARIMA aka ARFIMA(p,d,q) models. http://CRAN.R project.org/package=fracdiff. Fraser S., Dickson B.L., 2007. A new method for data integration and integrated data interpretation: self organizing maps. In Proceedings of Exploration 07: Fifth Decennial International Conference on Mineral Exploration, edited by B. Milkereit, 907 910. Friedel, M.J., 2011a. Climate change effects on ecosystem services in the United States Issues of national and global security. In: Baba, A., Tayfur, G., Howard, K.W.F., Friedel, M.J., Chambel, A., 2011, Climate Change and its Effect on Water Resources Issues of National and Global Security, NATO Science for Peace and Security Science Series Environmental Security, vol. 3, Springer, Dordrecht, The Netherlands, 318 p. ISBN 978 94 007 1142 6. Friedel, M.J., 2011b. A data driven approach for modeling post fire debris flow volumes and their uncert ainty, Environmental Modelling & Software 26 (12), 1583 1598. Friedel, M.J., 2012a, Data driven modeling of surface temperature anomaly and solar activity trends, Environmental Modelling & Software 37, 217 232. Friedel, M.J., 2012b, Hybrid modeling to pre dict the economic feasibility of mining undiscovered porphyry copper deposits. Applied Soft Computing 13, 1016 1032. Friedel., M.J., Asch, T., Oden, C ., 2012. Hybrid analysis of multiaxis electromagnetic data for discrimination of munitions and explosives of concern, Geophysical Journal International 190 (2), 960 980.
PAGE 58
50 Friedel, M. J. Akbari Esfahani, A., 2013 Influence of global volcanic processes on climate modulation in the southwestern US T24 presented at 201 3 Fall Meeting, GSA Denver C olorado ., 28 Oct. 1 Nov Furnis, M.J., 2010. Water, climate change, and forests. U.S. Forest service, General Technical Report, PNW GTR 812, pp. 75. Gray, L.J., Beer, J., Geller, M., Haigh, J.D., Lockwood, M., Matthes, K. Cubasch, U., Fleitmann, D., Harrison, G., Ho od, L., Luterbacher, J., Meehl, G.A., Shindell, D., van Geel, B., White, W., 2010. Solar influences on climate. Reviews of Geophysics 48: RG4001, 1 53. Haslett, J., Raftery, A.E., (1989). Space time Modelling with Long memory Dependence: Assessing Ireland' s Wind Power Resource (with Discussion); Applied Statistics 38, 1 50. Hastie, T., Tibshirani, R., Friedman, J., 2002. The elements of statistical learning, Springer Verlag, Berlin, 533 pp. Hornik, R., 2011. The R FAQ, ISBN 3 900051 08 9. http://www.r proje ct.org/, last visited 1 19 2012. Hurst, H.E., 1951. Long term storage capacity of reservoirs, Trans. Am. Soc. Civ. Eng., 116, 770 779, 1951 Hyndman, R.J., 2011. Forecast: Forecasting functions for time series. http://CRAN.R project.org/package=forecast. Hyndman, R.J., Koehler, A.B., 2006. Another look at measures of forecast accuracy, International Journal of Forecasting 22 (4), 679 688.
PAGE 59
51 Jain, A., Kumar, A.M., 2007. Hybrid neural network models for hydrologic time series forecasting. Applied Soft Computin g 7 (2), 585 592. Kalteh, A.M., Berndtsson, R., 2008. Review of the self organizing map (SOM) approach in water resources: Analysis, modeling and application. Environmental Modelling & Software 23, 835 845. Kaplan, D.T., 1994. Exceptional events as evidenc e for determinism. Physica D 73, 38 48. Koenker, R., 2005. Quantile regression, Econometric Society Monographs, Cambridge University Press, New York, 349 p. Kohonen T., 1984. Self Organization and Associative Memory. Springer, Berlin. Kohonen, T., 2001. Self organizing maps, Third Extended Edition, Springer Series in Information Sciences, Vol. 30, Springer, Berlin, Heidelberg, New York, 253 p. Loke, E., Arnbjerg box modelling: a comparison. In: Joliffe, I.B., Ball, J.E. (Eds.), Eighth International Conference: Urban Storm Drainage Proceedings, vol. 1. The Institution of Engineers Australia, Australia. Malek, M.A., Harun, S., Shamsuddin, S.M., Mohamad I., 2008. Imputation of time s eries data via Kohonen self organizing maps in the presence of missing data. Engineering and Technology 41, 501 506. Mandelbrot, B.B., 1967. How Long is the coast of Britain? statistical self similarity and fractional dimension, Science 156: 636 638, 1967. Mandelbrot, B.B., van Ness, J.W., 1968. Fractional brownian motions, fractional noises and applications, Siam Review 10 (4), 422 437.
PAGE 60
52 Mann, M.E., 2002. Little Ice Age. In Volume 1, The Earth system: Physical and chemical dimensions of global environmenta l change, Encyclopedia of Global Environmental Change (pp. 504 509). Chichester: John Wiley & Sons. Mann, M. E., Zhang, Z., Hughes, M. K., Bradley, R. S., Miller, S. K., Rutherford, S., Ni, F., 2008. Proxy based reconstructions of hemispheric and global su rface temperature variations over the past two millennia. Proceedings of the National Academy of Sciences, 105 (36), 13252 13257. Mann, M. E., Zhang, Z., Rutherford, S., Bradley, R. S., Hughes, M. K., Shindell, D., Ni, F., 2009. Global signatures and dynam ical origins of the little ice age and medieval climate anomaly. Science 326: 1256 1260, 2009. Martinez de Pison Ascacibar, F.J., Ordieres Meres, Joaquin, Castejon Limas, Manuel, de Cos Juez, F.J., 2008, Functions for calculating fractal dimension Package http://cran.r project.org/src/contrib/Archive/fdim/ National Oceanic and Atmospheric Administration, 2013. National Climatic Data Center. http://vlb.ncdc.noaa.gov/temp and precip/time series/ Nunes, S.A., Romani, L.A., Avila, A.M., Traina, .J., d e Sousa, E.P., Traina, A.J., 2011. Fractal based analysis to identify trend changes in multiple climate time series. Journal of Information and Data Management 2 (1), 51 57. Palmer, W.S., 1965. U.S. Weather Bureau Res. Pap. 45, 58 pp. Rae, W., Oxley, L., R eale, M., Brown, J., 2011. The emipirical properties of some popular estimators of long memory processes, in 19th International Congress on Modelling and Simulation, Perth, Australia, 2011.
PAGE 61
53 Rallo, R., Ferre Gine, J., Arenas, A., Giralt, F., 2002. Neural vi rtual sensor for the inferential prediction of product quality form process variables. Computers and Chemical Engineering 26 (12), 1735 54. Rampino, M.R., Self, S., 1982. Historic eruptions of Tambora (1815), Krakatau (1883), and Agung (1963), their strato spheric aerosols, and climatic impact. Quaternary Res. 18, 127 143. doi:10.1016/0033 5894(82)90065 5. Randhir, T., Ekness, P., 2009. Urbanization on watershed habitat potential: A multivariate assessment of thresholds and interactions. Ecohydrology 2, 88 1 01. Riebsame, W. E., Changnon, S.A., Karl, T.R., 1991. Drought and natural resources management in the United States: Impacts and implications of the 1987 89 drought. Westview Press New York, New York, 174 pp. Ritter, H. Schulten, K., 1986. On the station ary state of Kohonen's self organizing sensory mapping. Biological Cybernetics 54 (1), 99 106. Said, S., Dickey, D.A., 1984. Testing for unit roots in autoregressive moving average models of unknown order. Biometrika 71, 599 607. Schwing, F.B., Murphree, T., Green, P.M., 2002. The Northern Oscillation Index (NOI): A new climate index for the northeast Pacific. Progress in Oceanography 53, 115 139. http://www.esrl.noaa.gov/psd/gcos_wgsp/Timeseries/NAO/ Smith, T.M., Reynolds, R.W., 2003. Extended reconstruct ion of global sea surface temperatures based on COADS data (1854 1997). Journal of Climate 16, 1495 1510. Turcotte, D.L., 1997. Fractals and chaos in geology and geophysics, New York: Cambridge University Press.
PAGE 62
54 Ubilava, D. Helmers, .G., 2013. Forecastin g ENSO with a smooth transition autoregressive model, Environmental Modelling & Software 40, 181 190. ISSN 1364 8152, 10.1016/j.envsoft.2012.09.008. Vesanto J., Alhoniemi E., 2000. Clustering of the Self Organized Map. IEEE Transactions on Neural Networks 11(3), 586 600. Wang, S., 2003. Application of self organizing maps for data mining with incomplete data sets. Neural Computing and Applications 12, 42 48. Woodhouse, .A., Overpeck, J.T., 1998. 2000 years of drought variability in the central United State s. Bulletin of the American Meteorological Society 79 (12), 2693 2714.
