Marine Wind and Wave Height Trends at Different ERA-Interim Forecast Ranges
OLE JOHAN AARNES
Norwegian Meteorological Institute, and Geophysical Institute, University of Bergen, Bergen, Norway
SALEH ABDALLA,JEAN-RAYMOND BIDLOT, AND ØYVIND BREIVIK
European Centre for Medium-Range Weather Forecasts, Reading, United Kingdom
(Manuscript received 2 July 2014, in final form 16 September 2014)
ABSTRACT
Trends in marine wind speed and significant wave height are investigated using the global reanalysis ERA-
Interim over the period 1979–2012, based on monthly-mean and monthly-maximum data. Besides the tra-
ditional reanalysis, the authors include trends obtained at different forecast range, available up to 10 days
ahead. Any model biases that are corrected differently over time are likely to introduce spurious trends of
variable magnitude. However, at increased forecast range the model tends to relax, being less affected by
assimilation. Still, there is a trade-off between removing the impact of data assimilation at longer forecast
range and getting a lower level of uncertainty in the predictions at shorter forecast range. Because of the sheer
amount of assimilations made in ERA-Interim, directly and indirectly affecting the data, it is difficult, if not
impossible, to distinguish effects imposed by all updates. Here, special emphasis is put on the introduction of
wave altimeter data in August 1991, the only type of data directly affecting the wave field. From this, it is
shown that areas of higher model bias introduce quite different trends depending on forecast range, most
apparent in the North Atlantic and eastern tropical Pacific. Results are compared with 23 in situ measure-
ments, Envisat altimeter winds, and two stand-alone ECMWF operational wave model (EC-WAM) runs with
and without wave altimeter assimilation. Here, the 48-h forecast is suggested to be a better candidate for trend
estimates of wave height, mainly due to the step change imposed by altimeter observations. Even though wind
speed seems less affected by undesirable step changes, the authors believe that the 24–48-h forecast more
effectively filters out any unwanted effects.
1. Introduction
Long-term observation records of marine wind and
sea state are scarce in comparison with the time series
found over land routinely collected since the nineteenth
century in certain populated regions of Earth (Hurrell
1995). The need for long, reliable time series of marine
near-surface winds U
10
and significant wave height H
s
is
increasing as climate projections require a baseline cli-
matology against which to be compared, and even more
so if dynamical models of the sea state are to be included
in future coupled climate scenarios (Cavaleri et al. 2012;
Dobrynin et al. 2012; Sterl et al. 2012; Hemer et al. 2013;
Khon et al. 2014). There are also more immediate needs
for reliable time series of past wind and wave climate,
such as estimates of return values in areas without ob-
servational records (Caires and Sterl 2005; Aarnes et al.
2012; Breivik et al. 2013, 2014) or decadal trends in wind
and wave parameters. A number of recent regional studies
on wave climate variability and trends from hindcasts and
reanalyses are presented in Appendini et al. (2014)
(Gulf of Mexico); Reguero et al. (2013) and Izaguirre
et al. (2013) (Central and South America); Bromirski
et al. (2013) (North Pacific); Dodet et al. (2010),
Bertin et al. (2013), Wang and Swail (2001), and Wang
et al. (2012) (North Atlantic); and Wang and Swail
(2002) and Semedo et al. (2011) (Northern Hemi-
sphere). Most of these studies are to some extent di-
rectly or indirectly affected by data assimilation.
Therefore, trend estimates from altimeter wave heights
are an attractive alternative; see, for example, Hemer
et al. (2010) (Southern Hemisphere) and Young et al.
(2011, 2012) (global). These studies, however, rely
highly on intercalibration between satellite missions.
Corresponding author address: Ole Johan Aarnes, Norwegian
Meteorological Institute, Allégaten 70, NO-5007 Bergen, Norway.
E-mail: ole.aarnes@met.no
15 J
ANUARY 2015 A A R N E S E T A L . 819
DOI: 10.1175/JCLI-D-14-00470.1
2015 American Meteorological Society
Global marine wind and wave reanalyses can be pow-
erful proxies for observational records, provided that
they do indeed prove reliable in data-scarce regions and
periods. With the advent of long marine Earth-observing
satellite missions (European Remote Sensing Satellites
ERS-1 and ERS-2, starting in 1991), the quality and cov-
erage of marine U
10
and H
s
observations rose dramati-
cally. Since such observations are usually assimilated into
global reanalyses the question of the quality of a reanalysis
before and after the satellite era should be addressed.
Numerical weather prediction models are confined by
several characteristics, the inherent physics, resolution,
and numerics. Parameters such as U
10
and H
s
will vary
accordingly. If the modeled wind and wave climate
corresponds to the observed climate, assimilation be-
comes a mere way to keep the model on the right track.
Should the model however have an inherent bias rela-
tive to the observed climate, data assimilation will force
the model away from its own climate, only to revert back
as the model is integrated forward in time (forecast) less
constrained by assimilation. In this way, a model may end
up with different statistics at the time of analysis (here-
after ANA) and at increased forecast range (FCR). Not
surprisingly, ANA will correlate better with the real cli-
mate, the reason why reanalyses are particularly attractive
for s tudying past climate. However, as the number and
quality of assimilations vary with time, it is reasonable to
assume ANA to reflect these changes. Any discrepancies
found between the model output at ANA and increased
FCR may potentially identify any model biases, or
equivalently effects of assimilation. In the end, if we are
not correcting model bias at ANA, the trend should be
similar at ANA and increased FCR, assuming model drift
is insignificant.
ERA-Interim (hereafter ERA-I) is a reanalysis de-
veloped by the European Centre for Medium-Range
Weather Forecasts (ECMWF) under the European
Reanalysis Project (ERA), coupling atmosphere and
surface waves, covering the period from 1979 to the
present day (see Dee et al. 2011). Besides the traditional
reanalysis, ERA-I also performs 10-day forecasts from
reanalyzed fields twice a day. Since 1979, the number of
observations assimilated has increased substantially.
Between 1989 and 2010, there was a tenfold increase,
from 10
6
to 10
7
day
21
, with the biggest addition coming
from spaceborne instruments. The only wave data
products used in ERA-I are the significant wave height
observations from radar altimeters on board the satel-
lites ERS-1, ERS-2, Environmental Satellite (Envisat),
Jason-1, and Jason-2 available for different periods, but
with almost sustained existence since August 1991. Even
though the atmospheric model and the wave model are
coupled, the wave data analysis is done independently
from the atmospheric assimilation. The wave model
assimilation scheme, first developed by Lionello et al.
(1992), is based on the optimum interpolation (OI)
technique. This simple scheme provides an update to the
significant wave height field at ANA that is produced by
the coupled system following its four-dimensional vari-
ational assimilation of all sort of atmospheric data. This
update is then translated into an update of the model
spectra used for the subsequent model integration.
Ideally, a data assimilation scheme should only correct
for random error in the model; otherwise, the model
biases would quickly reappear once the model is in-
tegrated beyond ANA. This is even more prominent in
a forced system such as the wave model with an OI
scheme such as that employed in ERA-I.
Since altimeter wave height data originate from several
instruments with different characteristics and processing
procedures, there is a need for intercalibration to har-
monize them. For ERA-I, this was done with respect to in
situ observations (J. Bidlot 2015, unpublished manu-
script). Considering the fact that the wave-height climate
characteristics from in situ observations are not the same
as ERA-I model, first, ERA-I wave heights exhibit sys-
tematic biases with respect to in situ data (J. Bidlot 2015,
unpublished manuscript) and, consequently, biases with
respect to the intercalibrated altimeter data (Abdalla
et al. 2011). These biases are partially corrected at ANA
with the introduction of intercalibrated altimeter data.
However, before August 1991 or during periods when
altimeter data are temporarily not available, the ERA-I
wave heights return to their biased state.
In this study we investigate trends in H
s
and U
10
based
on ERA-I and how they are affected by nonstationary
assimilation over the period 1979–2012. The aim is to
detect any spurious trends and possibly propose an al-
ternative to ANA (i.e., an FCR offering trends more
representative of the real climate). Special attention is
being made on the era with and without altimeter wave
height assimilation. We compare ERA-I with observa-
tional records from buoys and satellites, taking care to
distinguish between independent datasets and datasets
that have been assimilated in the reanalysis. A similar
study based on atmospheric temperature from ERA-I
may be found in Simmons et al. (2014).
The paper is organized as follows. Section 2 introduces
the data [i.e., ERA-I, two stand-alone ECMWF opera-
tional wave model (EC-WAM) runs with and without
altimeter wave height assimilation] and observations, in
situ and altimetry. Section 3 presents the methodology for
trend estimation and the RHtestsV4 (Wang and Feng
2013) software package , a homogenization t ool used herein
to correct for step changes inherent in the in situ obser-
vations. Section 4 presents global trends in H
s
and U
10
and
820 JOURNAL OF CLIMATE VOLUME 28
comparisons between the different datasets. Section 5
discusses the results and relates the findings with similar
studies. Finally, section 6 offers some conclusions.
2. Data
a. ERA-Interim
ERA-I presents a third-generation reanalysis at the
ECMWF and possess a number of improvements from
its predecessors ERA-15 (Gibson et al. 1997) and ERA-
40 (Uppala et al. 2005); see Dee et al. (2011) for more
information. The ongoing project was originally meant
to improve the data-rich period of the 1990s and 2000s
following the appearance of Earth-observing satellites
such as ERS-1. In 2011 the reanalysis was extended
backwards from January 1989 to January 1979.
ERA-I is run with the same setup as the Integrated
Forecasting System (IFS) release cycle Cy31r2, used
operationally at ECMWF during the period December
2006 through June 2007. The horizontal resolution of the
atmospheric model is approximately 79 km (T255 spec-
tral truncation) on a reduced Gaussian grid. The cou-
pled wave component (Janssen 2004) is somewhat
coarser at approximately 110 km. The wave model is run
with shallow water physics where appropriate and dis-
cretized using 24 directions and 30 frequencies.
Prior to 2002, most of the observations assimilated in
ERA-I are similar to those used in ERA-40 with some
improvements [see section 4.1 of Dee et al. (2011) for
details]. In the context of this study, the most notable
changes were additional and recalibrated scatterometer
surface wind speeds, reprocessed wave data from ERS-1
and ERS-2 calibrated against buoy data and improved
satellite radiance data. After 2002, ERA-I uses obser-
vations from the ECMWF’s operational archive. ERA-I
assimilates altimeter H
s
data (see Table 1) and buoy U
10
data, while in situ H
s
and altimeter U
10
are not used in
the reanalysis and can serve to independently evaluate
the merits of the reanalyses.
ERA-I is run twice da ily, at 0000 and 120 0 UTC, but offers
6-hourly data at ANA, a blend of analysis and 6-hourly
reforecasts. Beyond the 12-h FCR, reforecasts are only
available at 0000 and 1200 UTC. To be consistent, we have
calculated monthly means and monthly maxima at ANA
and increased FCR based on 0000 and 1200 UTC only.
b. EC-WAM run with and without altimeter wave
height assimilation
In the following we compare two stand-alone EC-
WAM runs with and without wave altimeter assimila-
tion (hereafter denoted as WAM-AS and WAM-NAS,
respectively), that are otherwise identical. Both runs are
forced with archived ERA-I winds equivalent to U
ANA
10
and span the period 1992–2011 (20 yr). While ERA-I
uses a 30-min time step and evolving wind fields, the
wind fields forcing the stand-alone runs are sampled
every 6 h and kept constant over the same period. Un-
like ERA-I, which is run with two-way interaction, the
wave model is run in a separate operation with no
feedback from the wave model to the atmospheric
model. In addition, the stand-alone runs are based on
a later WAM cycle (Cy36r1). Since ERA-I (Cy31r2)
there have been three updates to the IFS WAM code
(see http://www.ecmwf.int/en/forecasts/documentation-
and-support/changes-ecmwf-model). In Cy33r1 (June
2008) the shallow water physics were improved by
modifying the nonlinear source term (Snl) and a new
advection scheme was implemented, reducing the gar-
den sprinkler effect (GSE). In Cy35r3 (September 2009)
the wave damping was intensified by including a weak
negative term in the wind input source term (Sin),
mainly affecting the longer wave components–reducing
swell. In Cy36r1 (January 2010) the wave model reso-
lution was increased from 0.368 to 0.258 . However, for
this particular experiment the resolution was set at 0.368,
going from ;110 to ;40 km. Further, the bathymetry
was refined from ETOPO5 (http://www.ngdc.noaa.gov/
mgg/global/etopo5.html) to ETOPO2 (http://www.ngdc.
noaa.gov/mgg/fliers/01mgg04.html), which better rep-
resents areas of partial blocking. The number of spectral
frequencies (30) and directions (24) are unchanged. The
WAM-AS was run with slightly improved bias correc-
tions of the altimeter data (i.e., ERA-I used to have
a nonoptimal bias correction up until the end of January
2010). Figure 1 presents the mean discrepancy in H
s
between ERA-I and WAM-AS over the period 1992–
2011, which illustrates the main features mentioned
TABLE 1. Altimeter H
s
assimilated in ERA-I.
ERS-1 1 Aug 1991–3 Jun 1996
ERS-2 3 May 1995–21 Jul 2003
Envisat 21 Jul 2003–April 2012
Jason-1 20 Oct 2003–3 Jul 2013
Jason-2 1 Feb 2010–present
FIG. 1. Mean discrepancy in H
s
between ERA-I (ANA) and
WAM-AS over the period 1992–2011.
15 J
ANUARY 2015 A A R N E S E T A L . 821
above. ERA-I is in general higher in the tropics because
of less swell damping and higher in the lee of poorly
resolved islands and shoals. The GSE is particularly
evident east-northeast of Hawaii.
c. Observations: In situ and Envisat
In an attempt to validate trends we use a selection of
23 in situ observations, H
s
and U
10
. We emphasize the H
s
data as they are not assimilated, but we also include U
10
for comparison. The hourly data are averaged over 62h
and centered around synoptic times (Bidlot et al. 2002).
Only collocated data with ERA-I are used to calculate
monthly means.
For the period November 2002–October 2010, a total
of nine years, we have binned Envisat altimeter winds
into 28328 latitude–longitude bins and collocated the
‘‘super observations’’ with U
ANA
10
in time and space. The
super observations represent altimeter data averaged
along track corresponding to the model resolution. Like
in situ H
s
, altimeter U
10
are not assimilated and there-
fore independent of ERA-I.
3. Method
a. Trend
As the probability density function (PDF) of H
s
and
U
10
do not conform to a Gaussian shape, trend estimates
should not be based on a simple regression. In the fol-
lowing, the magnitude of the trend is determined by
(Sen 1968; Yue et al. 2002)
Trend 5 median
x
j
2 x
l
j 2 l
" l , j, (1)
where x represents the data at time j and l. This offers
a robust estimate of any monotonic trend.
In the following, we apply the seasonal Kendall test,
a nonparametric test of randomness (H
0
) against trend
(H
1
), an extension of the Mann–Kendall test (Mann
1945; Kendall 1948) especially adapted to seasonal data
with serial dependence (Dietz and Killeen 1981; Hirsch
and Slack 1984). Let X
i
5 (x
i1
, x
i2
, ..., x
in
i
) represent
the monthly data of U
10
and H
s
, where n
i
is the total
number of entries from month i 5 1, 2, ..., 12 (only one
entry per year), then the seasonal Kendall statistics for
month i is expressed by
S
i
5
å
n
i
21
k51
å
n
i
j5k11
sgn(x
ij
2 x
ik
). (2)
In case of missing data at time j or k, sgn(x
ij
2 x
ik
) is set
to zero. From Mann (1945), Kendall (1948), and Hirsch
et al. (1982) we define S
0
5
å
12
i51
S
i
, having a mean and
variance given by
E[S
0
] 5
å
12
i51
E[S
i
] 5 0 and (3)
var[S
0
] 5
å
12
i51
var[S
i
] 1
å
12
i51
å
12
l51
cov(S
i
, S
l
) for i 6¼ l ,
(4)
where var[S
i
] 5 n
ig
(n
ig
2 1)(2n
ig
1 5)/18 and n
ig
represents
the number of nonmissing data per month (n
ig
5 n
i
for
complete series). According to Hirsch et al. (1982), cov
(S
i
, S
l
) 5 0 when S
i
and S
l
are independent random
variables. However, this fails to hold for monthly lag-1
serial correlation as low as 0.2 ( Hirsch and Slack 1984).
In the following we use an estimate of the covariance
term defined by Dietz and Killeen (1981), which is well
documented in Hirsch and Slack (1984). A two-sided
test for trend is based on the standard normal variate
Z defined by
Z 5
8
>
>
>
>
>
>
>
<
>
>
>
>
>
>
>
:
S
0
2 1
(var[S
0
])
1/2
,ifS
0
. 0
0, if S
0
5 0
S
0
1 1
(var[S
0
])
1/2
,ifS
0
, 0
, (5)
where H
0
is accepted when jZj , 1.96, using a signifi-
cance level of a 5 0.05.
As stated in Hirsch and Slack (1984) the covariance
term defined by Dietz and Killeen (1981) should be used
with some caution, such as for small sample sizes less
than 10 yr and in situations where data are, in fact, in-
dependent (e.g., in cases of many missing data). In the
following analysis we omit the covariance term when
found appropriate. This will be stated in the text.
b. RHtestsV4—Homogenization
The RHtestsV4 is a software package developed to
detect and adjust for sudden step changes, or shifts, in-
herent in time series for reasons other than climatic
changes. Here, we use the tool to homogenize monthly
in situ observations that have been altered because of,
for example, hardware and software updates. Even
though the RHtestsV4 is capable of detecting shifts by
analyzing the observed time series (referred to as base
series) solely by itself (Wang 2008a,b), Wang and Feng
(2013) highly recommend the use of a homogeneous and
well-correlated reference series for a more reliable re-
sult, especially when used in an automatic manner. In
the following, we use collocated ERA-I data. Shifts are
822 JOURNAL OF CLIMATE VOLUME 28
detected based on the penalized maximal t (PMT) test
(Wang et al. 2007; Wang 2008a), using the base-minus-
reference series, which is assumed to have zero-trend
and Gaussian errors.
Once any shifts have been established the different
segments of the base series are adjusted according to the
quantile-matching (QM) procedure presented in Wang
et al. (2010, their section 5) and Vincent et al. (2012).In
short, the empirical distributions of all detrended seg-
ments are matched. The software has shown promise in
Gemmrich et al. (2011) and Vincent et al. (2012).
It should be noted that, even though possible, we have
not run the RHtestsV4 with metadata describing buoy
updates, which may have increased the performance of
the RHtestsV4.
4. Results
We start off by investigating the mean difference in
monthly-mean H
s
and U
10
between ANA–24-h forecast
(ANA–FC24) and 24–48-h forecast (FC24–FC48) over
the periods 1979–91 and 1992–2012, that is, before and
after introducing wave height altimeter assimilation in
ERA-I (see Fig. 2). The red color scale in Fig. 2 shows
where H
s
and U
10
are decreasing with increased FCR,
while the blue color scale shows the opposite. The dis-
crepancies are mainly due to assimilation effects and
reflect model bias. For instance, in Fig. 2b, H
ANA
s
is
relatively higher in the storm track of the northeastern
Atlantic as the model has been corrected for a negative
bias related to insufficient wave growth in the area.
According to Hanley et al. (2010), the annual mean
wave age is lower in the North Atlantic storm track as
opposed to the North Pacific, related to a stronger mean
intensity in the extratropical cyclones (Bengtsson et al.
2006) and a corresponding stronger wind climate (Sterl
and Caires 2005). Further, waves are more fetch limited
in the North Atlantic because of a more varied coastline,
particularly near the southern tip of Greenland and in
the lee of Iceland. Fetch will also vary with season re-
sulting from ice extent. It is striking that the biggest
difference between H
ANA
s
and H
FC24
s
seems to follow the
mean ice edge of the winter season. In contrast, waves
generated in the storm tracks of the Southern Hemi-
sphere are far less fetch limited by land. And, since the
mean wind direction is westerly, with a slight northerly
component (see Hanley et al. 2010, their Fig. 4), ice
extent is probably not affecting the mean wave growth in
the same sense.
In the eastern tropical Pacific, an area more or less
completely dominated by swell (Semedo et al. 2011), the
assimilation effect is reversed. The model overestimates
the presence of swell, so H
ANA
s
is corrected down and
therefore relatively lower than H
FC24
s
. In the following,
our prime concern is how the model bias is dealt with
over time. Comparing Figs. 2a and 2b, the latter period is
clearly more influenced by assimilation. The corre-
sponding plots made with U
10
, Figs. 2c and 2d, do not
show the same geographical differences between the
two periods. However, there seems to be a more uniform
strengthening of U
ANA
10
in the latter period. For FC24–
FC48 (Figs. 2e–h ), the differences between the two pe-
riods are far less pronounced, indicating that the model
is relaxing toward its own climate, neglecting bias cor-
rections made at ANA.
To further illustrate the effect of assimilation, we plot
the discrepancies in monthly-mean H
s
and U
10
in-
tegrated over the Southern Hemisphere (SH; .208S),
the tropics (208S–208N), and the Northern Hemisphere
(NH; .208N) between ANA–FC24 and FC24–FC48, see
Figs. 3a–c and 3g–i. All ice-covered areas have been
removed from the H
s
and U
10
data. During the 1980s the
discrepancy between ANA–FC24 (in red) is fairly stable
for both parameters, with an exception found in the NH
where U
ANA
10
is steadily increasing relative to U
FC24
10
.In
August 1991 there is an abrupt jump in H
ANA
s
in all re-
gions caused by the sudden introduction of altimeter
observations. The effect is most pronounced north of
208N, but clearly visible south of 208S. In the tropics the
effect is negative and slightly less distinct. After August
1991 there is a steady relative increase in H
ANA
s
outside
the tropics up until approximately 2006. This seems re-
lated to a similar behavior in U
ANA
10
globally. Given the
fact that wave conditions in the tropics are highly
influenced by waves (swell) generated in the SH and
NH, it is somewhat surprising to find a relative decrease
in H
ANA
s
in the tropics, especially since the ‘‘local’’
U
ANA
10
also is increasing. In the last part of the period
H
ANA
s
and U
ANA
10
are decreasing. Overall, the largest
fluctuations in ANA–FC24 are found in U
10
in the
tropics, reflecting the poorer predictability in the area.
Notice that the FC24–FC48 comparison (Figs. 3g–i in
green) shows less fluctuation compared to ANA–FC24.
This is as expected as the model is gradually becoming
less influenced by assimilation with FCR. It should be
added that U
ANA
10
also shows evidence of minor step
changes. In 2000, when QuikSCAT was first introduced,
there is a clear drop in U
ANA
10
south of 208S (see Fig. 3g).
In Figs. 3d–f and 3j–l we present corresponding linear
trends in H
s
and U
10
as obtained at different FCR.
Again, trends are based on spatial averages over each
region per month and only data from ice-free areas have
been used. Globally, all trends are positive and signifi-
cant (not shown). The strongest trends are obtained at
ANA and decreasing with FCR. Trends in H
s
are rela-
tively stronger compared to U
10
. South of 208S trends
15 JANUARY 2015 A A R N E S E T A L . 823
