Predictability of extreme values in geophysical models

Extreme value theory in deterministic systems is concerned with unlikely large (or small) values of an observable evaluated along evolutions of the system. In this paper we study the finite-time predictability of extreme values, such as convection, energy, and wind speeds, in three geophysical models. We study whether finite-time Lyapunov exponents are larger or smaller for initial conditions leading to extremes. General statements on whether extreme values are better or less predictable are not possible: the predictability of extreme values depends on the observable, the attractor of the system, and the prediction lead time.


Introduction
Extreme value theory (EVT) for time series {X i } ∞ i=1 studies the limiting distribution of the partial maxima max(X 1 , . . . , X n ) as n → ∞. EVT was originally developed for time series of near-independent random variables (Leadbetter et al., 1983;Coles, 2001;Beirlant et al., 2004), but in the last decade EVT has been extended to chaotic deterministic dynamical systems (Haiman, 2003;Freitas and Freitas, 2008;Freitas et al., 2010;Gupta, 2010;Holland et al., 2012a,b;Felici et al., 2007a,b). In the latter context one considers time series generated by evaluating a scalar observable along evolutions of the system.
In this paper we study the predictability of extreme values, rather than their probability distribution. Lyapunov exponents measure predictability, but they are not very useful in forecasting applications. Firstly, Oseledec (1968) proved that Lyapunov exponents are almost everywhere constant functions of the initial condition. This means that all forecast outcomes are equally (un)predictable regardless of the current state of the system. Secondly, forecasts are typically made ahead for only a few days up to weeks, whereas Lyapunov exponents are asymptotic quantities, which are computed for time tending toward infinity.
Finite-time Lyapunov exponents (FTLEs) measure the exponential growth rate of nearby trajectories over a finite time, and typically they strongly depend on the initial condition. Nese (1989) and Abarbanel et al. (1991) showed that FTLEs for dynamical systems such as the Lorenz-63 model and the Hénon map can be very different along various parts of the attractor. Prasad and Ramaswamy (1999) and Datta and Ramaswamy (2003) demonstrated that distributions of FTLEs in intermittent dynamical systems are non-Gaussian, asymmetric, and often have fat tails. Similarly, Lai (2007) found universal distributions for FTLEs in systems having an attractor consisting of two distinct components. For studies on how FTLEs converge to their infinite-time counterparts see Bailey et al. (1997), Ziehmann et al. (1999Ziehmann et al. ( , 2000, and references therein. In this paper we study whether initial conditions leading to extreme values have larger or smaller FTLEs and how this depends on the prediction lead time and the attractor of the system. Section 2 introduces our concept of FTLEs based on norms that are tailored to the observable of interest, such as energy, convection, or wind speed. In Sect. 3 we study the distribution of FTLEs as a function of lead time and threshold for three geophysical models: a spectral truncation of the barotropic vorticity equation and the Lorenz-63 and Lorenz-96 models. Section 4 concludes the paper with a discussion.    Fig. 1: Distributions of the set Λ τ,q for the Lorenz-63 model and observable φ con : box plots as a function of lead time τ for different thresholds q. Boxes indicate the interquartile range, whiskers extend to the minimum and maximum values, and circles mark the median. For each lead time four thresholds are chosen: q = −∞ (gray) and the 95th, 97th, and 99th percentiles of the observable (resp. red, green, and blue).

Extreme events
The setting is a system of ordinary differential equationṡ and we denote the evolution operator by Φ t . In addition, we consider a scalar observable φ : R n → R. We define the extreme value domain as the set of points in phase space for which the observable exceeds a threshold q: We speak of an extreme event whenever a trajectory Φ t (x) enters the set E q . The papers on extreme value theory for dynamical systems cited in the introduction focused on the case where φ has the form φ(x) = g(dist(x, x)), where g : [0,∞) → R is a continuous function and x is a point on the attractor. However, as pointed out by Holland et al. (2012b), physically relevant observations are not of this form. For example, geophysical quantities such as vorticity or components of the velocity vector can often be written as where q ∈ R n is a nonzero vector. Positive quantities, such as energy, enstrophy, or squared wind speeds, can often be written as where Q = Q ⊤ is a nonnegative definite n× n matrix. In this paper we exclusively focus on observations of the form (2) and (3).

Finite-time Lyapunov exponents
Infinitesimal errors in the initial condition x of the system (1) in the direction v evolve as L t (x)v, where the matrix L t (x) is the solution of the matrix intial value probleṁ which is often referred to as the tangent linear equation or first variational equation. The singular vectors of L t (x) give Fig. 1. Distributions of the set τ,q for the Lorenz-63 model and observable φ con : box plots as a function of lead time τ for different thresholds q. Boxes indicate the interquartile range, whiskers extend to the minimum and maximum values, and circles mark the median. For each lead time, four thresholds are chosen: q = −∞ (grey) and the 95th, 97th, and 99th percentiles of the observable (red, green, and blue, respectively).

Methodology
In summary, our methodology consists of checking whether initial conditions leading to extreme values typically have larger or smaller finite-time Lyapunov exponents.

Extreme events
The setting is a system of ordinary differential equationṡ and we denote the evolution operator by t . In addition, we consider a scalar observable φ : R n → R. We define the extreme value domain as the set of points in phase space for which the observable exceeds a threshold q: We speak of an extreme event whenever a trajectory t (x) enters the set E q . The papers on extreme value theory for dynamical systems cited in the introduction focused on the case where φ has the form where g : [0, ∞) → R is a continuous function and x is a point on the attractor. However, as pointed out by Holland et al. (2012b), physically relevant observations are not of this form. For example, geophysical quantities such as vorticity or components of the velocity vector can often be written as where q ∈ R n is a nonzero vector. Positive quantities, such as energy, enstrophy, or squared wind speeds can often be written as where Q = Q is a nonnegative definite n × n matrix. In this paper we exclusively focus on observations of the form Eqs. (2) and (3).

Finite-time Lyapunov exponents
Infinitesimal errors in the initial condition x of the system Eq.
(1) in the direction v evolve as L t (x)v, where the matrix L t (x) is the solution of the matrix intial value probleṁ which is often referred to as the tangent linear equation or first variational equation. The singular vectors of L t (x) give the direction and magnitude of maximal error growth along a finite-time segment of an orbit in state space (Buizza and Palmer, 1995;Palmer et al., 1998;Leutbecher and Palmer, 2008).
To measure the growth of perturbations over a time interval [0, τ ], we use different norms at t = 0 and t = τ . Let C be the initial error covariance matrix, and let P be an orthogonal projection matrix (i.e. P = P 2 = P ). Let · denote the Euclidean norm. The Rayleigh-Ritz quotient is maximal if and only if v is the right singular vector of P L τ (x)C 1/2 associated with the largest singular value. We define the time-τ finite-time Lyapunov exponents (FTLE) as where σ 1 (x, τ ) ≥ σ 2 (x, τ ) ≥ · · · ≥ σ n (x, τ ) are the singular values of P L τ (x)C 1/2 . In the following we assume that the initial error covariance C is the identity matrix. We choose the projection matrix P such that the singular vectors will point in the direction of maximal growth of the observable φ. For an observable of the form Eq.
(2), we set In case of observable Eq. (3), recall that a symmetric nonnegative definite n × n matrix has k ≤ n positive eigenvalues µ 1 ≥ · · · ≥ µ k > 0 and the corresponding eigenvectors q i are orthonormal. We can rewrite Eq. (3) as In this case we set Hence, the φ(x) is large if and only if P x is large.

Predictability of extreme values
We fix an initial point x 0 on the attractor and a sampling frequency ω > 0 and compute a sample of points along an orbit on the attractor: We choose ω and N sufficiently large to provide a good sampling of the attractor as far as both local and global fluctuations are concerned. Typically, we take N = 10 6 , whereas ω depends on the time scale of the system of interest.
To study the predictability of extreme events, we define the set containing the FTLEs for initial conditions such that after τ units of time the observable exceeds the threshold q. We study how the distribution of this set changes as the lead time τ and threshold q increase. For fixed lead times we compare the distributions of the set τ,q for q = −∞ (the entire sample) and for q being the 95th, 97th, and 99th percentiles of the sample If FTLEs cluster to larger positive (smaller negative) values when q increases, then extremes are less (better) predictable.

Results
We apply the methodology of Sect. 2 to three geophysical models: the Lorenz-63 and Lorenz-96 models and a spectral truncation of the barotropic vorticity equation. As observations we will take energy, convection, and squared wind speeds.

The Lorenz-63 model
In this section we consider the classical model for Rayleigh-Bénard convection derived by Lorenz (1963): with the classical parameter values σ = 10, ρ = 28, and β = 8/3. We use the observables φ con (x, y, z) = x and φ en (x, y, z) = x 2 + y 2 + z 2 that respectively measure the intensity of convection and the total energy. For the projection matrix in Eq. (5), we respectively take  = xy − βz, (6) , ρ = 28, and β = x 2 + y 2 + z 2 f convection and in (5) we respec-1 0 0 0 1 0 0 0 1   mple of N = 10 6 equency ω = 100. bles evaluated at set Λ τ,q changes ervable φ con . For l conditions leade have a negative of convection are or τ > 0.25 more leading to an ex-E which is larger that extreme valhen the lead time fferent thresholds difference in premes dissappears. set Λ τ,q changes or the observable  Black dots mark initial conditions leading to an exceedance of the 95th percentile of the observable φ con . For τ = 0.25, initial conditions leading to extremes typically have low FTLE, but for τ = 0.25, they typically have a large FTLE. This difference vanishes for longer lead times.
By numerical integration we computed a sample of N = 10 6 points on the attractor using the sample frequency ω = 100. Table 1 lists some statistics of the observables evaluated at the sample points. Figure 1 shows how the distribution of the set τ,q changes with lead time and threshold q for the observable φ con . For τ ≤ 0.25 more than 75 percent of the inital conditions leading to an exceedance of the 99th percentile have a Only for lea tions of Λ τ, shown). Fig. 3 and φ en the FTL various lead which lead t stated lead t large FTLE where the tw each other. regions quic conditions l coincide wit whereas init not.

The ba
In this sectio (BVE) for at where ψ is t Coriolis forc erator. In ad negative FTLE. This suggest that extreme values of convection are well-predictable up to small lead times. For τ > 0.25 more than 75 percent of the initial conditions leading to an exceedance of the 99th quantile have a FTLE that is larger than the median of τ,−∞ , which indicates that extreme values for this lead time are less predictable. When the lead time τ increases, the distributions of τ,q for different thresholds q become similar; this suggests that the difference in predictability between extremes and non-extremes disappears. the direction and magnitude of maximal error growth along a finite-time segment of an orbit in state space (Buizza and Palmer, 1995;Palmer et al., 1998;Leutbecher and Palmer, 2008).
To measure the growth of perturbations over a time interval [0,τ ] we use different norms at t = 0 and t = τ . Let C be the initial error covariance matrix, and let P be an orthogonal projection matrix (i.e., P = P 2 = P ⊤ ). Let · denote the Euclidean norm. The Rayleigh-Ritz quotient is maximal if and only if v is the right singular vector of P L τ (x)C 1/2 associated with the largest singular value. We define the time-τ finite-time Lyapunov exponents (FTLE) as where σ 1 (x,τ ) ≥ σ 2 (x,τ ) ≥ ··· ≥ σ n (x,τ ) are the singular values of P L τ (x)C 1/2 . In the following we assume that the initial error covariance C is the identity matrix. We choose the projection matrix P such that the singular vectors will point in the direction of maximal growth of the observable φ. For an observable of the form (2) we set In case of observable (3) recall that a symmetric nonnegative definite n × n matrix has k ≤ n positive eigenvalues µ 1 ≥ ··· ≥ µ k > 0 and the corresponding eigenvectors q i are orthonormal. We can rewrite (3) as In this case we set Hence, the φ(x) is large if and only if P x is large.

Predictability of extreme values
We fix an initial point x 0 on the attractor and a sampling frequency ω > 0 and compute a sample of points along an orbit on the attractor: We choose ω and N sufficiently large to provide a good sampling of the attractor as far as both local and global fluctuations are concerned. Typically, we take N = 10 6 , whereas ω depends on the time scale of the system of interest.
To study the predictability of extreme events we define the set containing the FTLEs for initial conditions such that after τ units of time the observable exceeds the threshold q. We Fig. 4. As Fig. 1, but for the observable φ en . Figure 4 shows how the distribution of the set τ,q changes with increasing lead time and threshold for the observable φ en . For all lead times τ almost all initial conditions leading to an exceedance of the 95th percentile of the energy have a FTLE that is larger than the median of τ,−∞ . This means that in this case extreme values of the energy observable φ en are systematically less predictable than non-extreme values.
Only for lead times much longer than τ = 1.75, the distributions of τ,q for different quantiles q become similar (not shown). Figures 2 and 3 show how for the observables φ con and φ en the FTLEs vary on different parts on the attractor for various lead times. Black dots mark the initial conditions which lead to an exceedance of the 95th percentile within the stated lead time. For short lead times, initial conditions with large FTLE are located in the central region of the attractor where the two wings of the butterfly-shaped attractor cross each other. When the lead time increases the unpredictable regions quickly spread out on the entire attractor. The initial conditions leading to extreme energy values systematically coincide with the parts on the attractor having larger FTLEs, whereas initial conditions leading to extreme convection do not.

The barotropic vorticity equation
In this section we consider the barotropic vorticity equation (BVE) for atmospheric flow over orography: where ψ is the stream function, ψ * is the forcing, βy is the Coriolis force, h is the orography, and is the Laplace operator. In addition, J is the Jacobian operator: A are systematically less predictable than non-extreme values.
Only for lead times much longer than τ = 1.75 the distributions of Λ τ,q for different quantiles q become similar (not shown). Fig. 3 and Fig. 4 show how for the observables φ con and φ en the FTLEs vary on different parts on the attractor for various lead times. Black dots mark the initial conditions which lead to an exceedance of the 95th percentile within the stated lead time. For short lead times initial conditions with large FTLE are located in the central region of the attractor where the two wings of the butterfly-shaped attractor cross each other. When the lead time increases the unpredictable regions quickly spread out on the entire attractor. The inital conditions leading to extreme energy values systematically coincide with the parts on the attractor having larger FTLEs, whereas initial conditions leading to extreme convection do     Fig. 6: As Fig. 1, but for the model (8)  which respectively mean that on the lateral boundaries the northward component of the velocity vanishes and that the net zonal flow is zero. Charney and DeVore (1979) studied a 3-dimensional spectral truncation of (7), and detected coexisting equilibria representing blocked and zonal flow patterns. De Swart (1989) studied higher order truncations, including the 6-dimensional truncatioṅ We study this model with the parameters of Crommelin et al. (2004), which are listed in Table 2. With this choice the channel [0,2π]×[0,πb] has dimensions 5000×1250 km, the orography h(x,y) = cos(x)sin(y/b) has a dimensional amplitude of 200 m, and one unit of time corresponds to one day. Crommelin et al. (2004) found intermittent transitions between two regions in state space corresponding to zonal and blocked flows, see Fig. 5. The dynamics consists of three recurrent episodes: (1) transitions from zonal to blocked flows, typically taking 30 days, (2) transitions from blocked to zonal flows, typically taking 40-80 days, and (3) spiraling behaviour around the zonal regime, typically lasting more than 200 days. This intermittent behaviour can be explained by the presence of Shil'nikov-like strange attractors due to homoclinic bifurcations taking place near a Hopf-saddle-node bifurcation (Broer and Vegter, 1984;Broer and Vitolo, 2008).
As observable we use the squared wind speed at several locations in the β-plane channel. In Appendix A we show that the east-and northward components of the wind vector at (x,y) can be computed from the state vectorψ as q ⊤ 1ψ and q ⊤ 2ψ , where q 1 ,q 2 ∈ R 6 are orthogonal vectors depending on (x,y). Hence, the squared wind speed at location (x,y) is given by which respectively mean that on the lateral boundaries the northward component of the velocity vanishes and that the net zonal flow is zero. Charney and DeVore (1979) studied a 3-dimensional spectral truncation of Eq. (7), and detected coexisting equilibria representing blocked and zonal flow patterns. De Swart (1989) studied higher order truncations, including the 6dimensional truncatioṅ We study this model with the parameters of Crommelin et al. (2004), which are listed in Table 2. With this choice the channel [0, 2π] × [0, πb] has dimensions 5000 × 1250 km, the orography h(x, y) = cos(x) sin(y/b) has a dimensional amplitude of 200 m, and one unit of time corresponds to one day. Crommelin et al. (2004) found intermittent transitions between two regions in state space corresponding to zonal and blocked flows, see Fig. 5. The dynamics consists of three recurrent episodes: (1) transitions from zonal to blocked flows, typically taking 30 days, (2) transitions from blocked to zonal flows, typically taking 40-80 days, and (3) spiraling behaviour around the zonal regime, typically lasting more than 200 days. This intermittent behaviour can be explained by the presence of Shil'nikov-like strange attractors, due to homoclinic bifurcations taking place near a Hopf-saddle-node bifurcation (Broer and Vegter, 1984;Broer and Vitolo, 2008).
As observable we use the squared wind speed at several locations in the β-plane channel. In Appendix A we show that the east-and northward components of the wind vector at (x, y) can be computed from the state vectorψ as q 1ψ and q 2ψ , where q 1 , q 2 ∈ R 6 are orthogonal vectors depending on (x, y). Hence, the squared wind speed at location (x, y) is given by φ(ψ) = (q 1ψ ) 2 + (q 2ψ ) 2 =ψ (q 1 q 1 + q 2 q 2 )ψ. For the projection matrix in Eq. (5) we take P = q 1 q 1 q 1 2 + q 2 q 2 q 2 2 .
We will study the extreme wind speeds at three locations: the side of the orography (x, y) = (0, 3π b/4), the top (x, y) = (0, π b/2), and the valley (x, y) = (π, π b/2), see Fig. 5. The corresponding wind speed observables are denoted by φ side , φ top , and φ val . By numerical integration we computed a sample of N = 10 6 points on the attractor using the sample frequency ω = 2. Table 1  For the projection matrix in (5) we take We will study the extreme wind speeds at three locations: the side of the orography (x,y) = (0,3πb/4), the top (x,y) = (0,πb/2), and the valley (x,y) = (π,πb/2), see Fig. 5. The corresponding wind speed observables are denoted by φ side , φ top , and φ val . By numerical integration we computed a sample of N = 10 6 points on the attractor using the sample frequency ω = 2. Table 1 lists some statistics of the observables evaluated at the sample points. Fig. 6, 7, and 8 show how the distributions of the set Λ τ,q change with lead time τ and threshold q for three observables φ sid , φ top , and φ val . For φ sid the initial conditions leading to extremes systematically have a large FTLE, which means that extremes are less predictable. For φ top the FTLEs for initial conditions leading to extremes are highly variable. In this case one cannot say that extremes are either better or less predictable than generic events. Finally, for φ val the FTLEs for initial conditions leading to extremes are systematically low: more than 75% of the initial conditions leading to extremes have an FTLE smaller than the median FTLE for all events.  For the projection matrix in (5) we take We will study the extreme wind speeds at three locations: the side of the orography (x,y) = (0,3πb/4), the top (x,y) = (0,πb/2), and the valley (x,y) = (π,πb/2), see Fig. 5. The corresponding wind speed observables are denoted by φ side , φ top , and φ val . By numerical integration we computed a sample of N = 10 6 points on the attractor using the sample frequency ω = 2. Table 1 lists some statistics of the observables evaluated at the sample points. Fig. 6, 7, and 8 show how the distributions of the set Λ τ,q change with lead time τ and threshold q for three observables φ sid , φ top , and φ val . For φ sid the initial conditions leading to extremes systematically have a large FTLE, which means that extremes are less predictable. For φ top the FTLEs for initial conditions leading to extremes are highly variable. In this case one cannot say that extremes are either better or less predictable than generic events. Finally, for φ val the FTLEs for initial conditions leading to extremes are systematically low: more than 75% of the initial conditions leading to extremes have an FTLE smaller than the median FTLE for all events.  Figures 6, 7, and 8 show how the distributions of the set τ,q change with lead time τ and threshold q for three observables φ sid , φ top , and φ val . For φ sid the initial conditions leading to extremes systematically have a large FTLE, which means that extremes are less predictable. For φ top the FTLEs for initial conditions leading to extremes are highly variable. In this case one cannot say that extremes are either better or less predictable than generic events. Finally, for φ val the FTLEs for initial conditions leading to extremes are systematically low: more than 75% of the initial conditions leading to extremes have an FTLE smaller than the median FTLE for all events. Figures 9-11 show the FTLEs vary on different parts of the attractor. For all observables the initial conditions with large FTLE are located near the zonal regime and the initial conditions with small FTLE are located near the blocked regime. This separation persists remarkably for increasing lead time. This is very different from the Lorenz-63 attractor, for which the unpredictable regions quickly spread out along the attractor. For the observables φ sid and φ top the initial conditions leading to extreme wind speeds systematically coincide with the unpredictable parts of the attractor. For φ val the initial conditions leading to extremes intersect the more predictable regions on the attractor.

The Lorenz-96 model
The Lorenz-96 model is a conceptual model for traveling waves in the atmosphere (Lorenz, 1996). This model is used to test data assimilation methods (Trevisan and Palatella, 2011) and subgrid scale parameterizations (Crommelin and Vanden-Eijnden, 2008). Although the model is not derived from physical principles it has features which are typical for geophysical models: forcing, dissipation, and energy preserving quadratic terms.
The model variables x 1 , . . . , x n can be interpreted as meteorological quantities, such as pressure or vorticity, along a circle of constant latitude where the index i of each variable x i plays the role of longitude. The dynamical equations are This separation persists remarkably for increasing lead time. This is very different from the Lorenz-63 attractor, for which the unpredictable regions quickly spread out along the attractor. For the observables φ sid and φ top the initial conditions leading to extreme wind speeds systematically coincide with the unpredictable parts of the attractor. For φ val the initial conditions leading to extremes intersect the more predictable regions on the attractor.

The Lorenz-96 model
The Lorenz-96 model is a conceptual model for traveling waves in the atmosphere (Lorenz, 1996). This model is used to test data assimilation methods (Trevisan and Palatella, 2011) and subgrid scale parameterizations (Crommelin and Vanden-Eijnden, 2008). Although the model is not derived from physical principles it has features which are typical for geophysical models: forcing, dissipation, and energy preserving quadratic terms. The model variables x 1 ,...,x n can be interpreted as meteorological quantities, such as pressure or vorticity, along a circle of constant latitude, where the index i of each variable x i plays the role of longitude. The dynamical equations are with periodic 'boundary conditions' x i+n = x i . The dimension n and forcing F > 0 are free parameters. We set n = 36 and F = 8 which were orginally used by Lorenz.
We use the observable x 2 i Fig. 9. Projection on the ( ψ 3 , ψ 1 )-plane of a the attractor of Eq. (8) with colour indicating the magnitude of the maximal time-τ FTLE. The grey circle (square) indicates the point in phase space corresponding with the zonal (blocked) regime of Fig. 5. Black dots mark initial conditions leading to an exceedance of the 95th percentile of the observable φ sid . Note that initial conditions leading to high wind speeds are located near the zonal regime, where also the FTLEs are large.
with periodic "boundary conditions" x i+n = x i . The dimension n and forcing F > 0 are free parameters. We set n = 36 and F = 8, which were originally used by Lorenz.
We use the observable x 2 i which measures the total energy of the system. For the projection in Eq. (5), we take the identity matrix. By numerical integration we computed a sample of N = 10 6 points on the attractor using the sample frequency ω = 100. Table 1   8 A This separation persists remarkably for increasing lead time. This is very different from the Lorenz-63 attractor, for which the unpredictable regions quickly spread out along the attractor. For the observables φ sid and φ top the initial conditions leading to extreme wind speeds systematically coincide with the unpredictable parts of the attractor. For φ val the initial conditions leading to extremes intersect the more predictable regions on the attractor.

The Lorenz-96 model
The Lorenz-96 model is a conceptual model for traveling waves in the atmosphere (Lorenz, 1996). This model is used to test data assimilation methods (Trevisan and Palatella, 2011) and subgrid scale parameterizations (Crommelin and Vanden-Eijnden, 2008). Although the model is not derived from physical principles it has features which are typical for geophysical models: forcing, dissipation, and energy preserving quadratic terms. The model variables x 1 ,...,x n can be interpreted as meteorological quantities, such as pressure or vorticity, along a circle of constant latitude, where the index i of each variable x i plays the role of longitude. The dynamical equations are with periodic 'boundary conditions' x i+n = x i . The dimension n and forcing F > 0 are free parameters. We set n = 36 and F = 8 which were orginally used by Lorenz.
We use the observable x 2 i Fig. 10. As Fig. 9, but for observable φ top . In this case initial conditions leading to large wind speeds are not confined to the unpredictable region of the attractor.
lists some statistics of the observable evaluated at the sample points. Figure 12 shows how the distribution of the set τ,q changes with lead time τ and threshold q. For all lead times τ , the distribution of the maximal FTLE does not change very much with increasing threshold q. This suggests that extreme values of the energy observable are neither better nor less predictable than non-extreme values.

Conclusions and discussion
In this paper we investigated the predictability of extreme values in geophysical models. We studied how (b) τ = 2.5 Fig. 11: As Fig. 9, but for observable φ val . In this case initial conditions leading to large wind speeds are located in the more predictable regions of the attractor.
which measures the total energy of the system. For the projection in (5) we take the identity matrix. By numerical integration we computed a sample of N = 10 6 points on the attractor using the sample frequency ω = 100. Table 1 lists some statistics of the observable evaluated at the sample points. Fig. 12 show how the distribution of the set Λ τ,q changes with lead time τ and threshold q. For all lead times τ the distribution of the maximal FTLE does not change very much with increasing threshold q. This suggests that extreme values of the energy observable are neither better nor less predictable than non-extreme values.

Conclusions and discussion
In this paper we investigated the predictability of extreme values in geophysical models. We studied how FTLEs depend on forecast lead time and the threshold on the observ-

Model
Obs. Predictability of extremes Lorenz-63 φcon For lead times up to τ = 0.25 extremes are well-predictable (negative FTLEs), but for longer lead times they are just as predictable as non-extremes.
φen For lead times up to τ = 01.75 extremes are less predictable than non-extremes. BVE φ sid Extremes are less predictable, and this persists up to long lead times due to the intermittent nature of the dynamics.
φtop Extremes are better predictable than nonextremes.
φ val Extremes are neither better nor less predictable than non-extremes.  Extremes are neither better nor less predictable than non-extremes. able. General statements on whether extreme values are better or worse predictable are not possible. Whether initial conditions leading to extreme values have larger FTLEs depends on (1) the observable, (2) the attractor of the system, and (3) the prediction lead time. Table 3 presents the main conclusion for each model/observable pair. FTLEs measure forecast error growth rate under the assumption that errors in the initial condition are infinitesimally small. For finite-size, but small, errors the maximal FTLE still is a good estimate of error growth (Harle et al., 2006). For larger errors, however, error growth may no longer be exponential and this will also effect the predictability of extreme values. The present study should therefore be extended to the setting of ensemble forecasts in which the predictability of extremes is measured in terms of the dispersion of ensemble members. This approach would also be more appropriate for operational weather forecasting models without a tangent linear model needed to compute FTLEs, such as the MOGREPS system of the UK Met Office (Bowler et al., 2008).
Another important question is: how predictable are realworld extremes, such as wind storms? Large-scale flow patterns, such as the North Atlantic Oscillation, cause temporal clustering of storms (Mailier et al., 2006;Vitolo et al., 2009). The emergence of these patterns might be a manifestation of intermittency, i.e., the irregular alternation between phases of chaotic and non-chaotic, such as steady or periodic, dynamics (Pomeau and Manneville, 1980). For example, the spectral truncation (8) of the barotropic vorticity equation exhibits intermittent transitions between zonal and blocked flows due to Shil'nikov-like strange attractors appearing near a Hopf-saddle-node bifurcation (Crommelin et al., 2004;Broer and Vitolo, 2008;Broer and Vegter, 1984). Other forms of intermittency due to bifurcations of planetary Fig. 11. As Fig. 9, but for observable φ val . In this case initial conditions leading to large wind speeds are located in the more predictable regions of the attractor. conditions leading to extreme values have larger FTLEs depends on (1) the observable, (2) the attractor of the system, and (3) the prediction lead time. Table 3 presents the main conclusion for each model/observable pair.
FTLEs measure forecast error growth rate under the assumption that errors in the initial condition are infinitesimally small. For finite-size, but small, errors, the maximal FTLE still is a good estimate of error growth (Harle et al., 2006). For larger errors, however, error growth may no longer be exponential and this will also effect the predictability of extreme values. The present study should, therefore, be extended to the setting of ensemble forecasts in which the predictability of extremes is measured in terms of the dispersion of ensemble members. This approach would also be more appropriate for operational weather forecasting models without a tangent linear model needed to compute FTLEs, such as the MOGREPS system of the UK Met Office (Bowler et al., 2008). Another important question is: how predictable are realworld extremes, such as wind storms? Large-scale flow patterns, such as the North Atlantic Oscillation, cause temporal clustering of storms (Mailier et al., 2006;Vitolo et al., 2009). The emergence of these patterns might be a manifestation of intermittency, i.e. the irregular alternation between phases of chaotic and non-chaotic, such as steady or periodic dynamics (Pomeau and Manneville, 1980). For example, the spectral truncation Eq. (8) of the barotropic vorticity equation exhibits intermittent transitions between zonal and blocked flows due to Shil'nikov-like strange attractors appearing near a Hopf-saddle-node bifurcation (Crommelin et al., 2004;Broer and Vitolo, 2008;Broer and Vegter, 1984). Other forms of intermittency due to bifurcations of planetary waves have been detected in low-order models of the shallow water equations (Sterk et al., 2010). Because different intermittent phases can have different error growth rates, the emergence of large-scale flow patterns might enhance predictability of extremes.
Finally, we note that FTLEs might not be the best measures of finite-time predictability. A different approach would be to apply techniques based on Takens' reconstruction theorem, such as correlation integrals and entropy; see Broer and Takens (2011) and references therein. Such techniques have been applied to develop early warning systems for thermal excursions in chemical reactors (Zaldívar et al., 2005). Potentially, these techniques can be useful in the prediction of extreme values in geophysical applications. We believe that these questions and problems will have sufficient potential for future research. waves have been detected in low-order models of the shallow water equations (Sterk et al., 2010). Because different intermittent phases can have different error growth rates, the emergence of large-scale flow patterns might enhance predictability of extremes.
Finally, we note that FTLEs might not be the best measures of finite-time predictability. A different approach would be to apply techniques based on Takens' reconstruction theorem, such as correlation integrals and entropy; see Broer and Takens (2011) and references therein. Such techniques have been applied to develop early warning systems for thermal excursions in chemical reactors (Zaldívar et al., 2005). Potentially, these techniques can be useful in the prediction of extreme values in geophysical applications. We believe that these questions and problems will have sufficient potential for future research.
De Swart (1988)  Note that the vectors q 1 and q 2 are orthogonal.
Acknowledgements. This research is part of the project "PREDEX: Predictability of Extreme Weather Events", funded by Complexity-NET: (www.complexitynet.eu). The authors gratefully acknowledge support by the UK and Dutch funding agencies involved in Fig. 12. As Fig. 1, but for the Lorenz-96 model Eq. (9) with observable φ en .
De Swart (1988) uses the convention u = −ψ y , v = ψ x . This gives the following expressions for the east-and northward components of the velocity field: u(x, y) = Note that the vectors q 1 and q 2 are orthogonal.