## Abstract

A nonlinear dynamical perspective on climate prediction is outlined, based on a treatment of climate as the attractor of a nonlinear dynamical system *D* with distinct quasi-stationary regimes. The main application is toward anthropogenic climate change, considered as the response of *D* to a small-amplitude imposed forcing **f**.

The primary features of this perspective can be summarized as follows. First, the response to **f** will be manifest primarily in terms of changes to the residence frequency associated with the quasi-stationary regimes. Second, the geographical structures of these regimes will be relatively insensitive to **f**. Third, the large-scale signal will be most strongly influenced by **f** in rather localized regions of space and time. In this perspective, the signal arising from **f** will be strongly dependent of *D*’s natural variability.

A theoretical framework for the perspective is developed based on a singular vector decomposition of *D*’s tangent propagator. Evidence for the dyamical perspective is drawn from a number of observational and modeling studies of intraseasonal, interannual, and interdecadal variability, and from climate change integrations. It is claimed that the dynamical perspective might resolve the apparent discrepancy in global warming trends deduced from surface and free troposphere temperature measurements.

A number of specific recommendations for the evaluation of climate models are put forward, based on the ideas developed in this paper.

## 1. Introduction

This paper is an extension and update of an informal article (Palmer 1993) that discussed a nonlinear dynamical perspective on climate prediction. Although many of the ideas developed in the current paper are relevant to naturally occurring climate fluctuations, the main application is to anthropogenic climate change.

Of course, one might argue that conceptual perspectives on climate change are irrelevant, since general circulation models (GCMs) provide explicit quantitative estimates of such change. The problem with such an argument is that GCMs are not perfect replications of the real world, and thus GCM predictions of climate change are themselves subject to errors and uncertainties in model formulation. The extent to which a model can be imperfect, yet useful for forecasting climate change, does depend on one’s conceptual perspectives about climate, particularly the linearity of climate. For example, a simplistic linear perspective would hold that climate change forecasts could be corrected for model imperfections by taking differences between a control integration (with present day CO_{2}), and a perturbed integration (with projected estimate of CO_{2}). A related perspective might suggest that one can tolerate errors in simulating modes of interannual variability, provided one is only interested in climate change on decadal and longer timescales. Here linearity is invoked to justify some sort of temporal filtering. Yet another linear perspective is that spatially localized errors in GCM performance can be tolerated provided one is only interested in the impact of anthropogenic forcing on ultra-large scales, such as hemispheric-mean temperature.

There are other examples where the underlying paradigm is of importance. The value of flux correction techniques depends on the extent to which misrepresentation of the nonlinear equations of motion can be mitigated by the addition of empirical corrective terms that force the model to a prescribed state (e.g., associated with present-day mean climate). In reality, this prescribed state is not itself a time-invariant solution of the equations. However, in the model it is possible that these corrective terms may excessively stabilize the prescribed state, leading to a damping of (low-frequency) variability about this state. In a nonlinear world, the transient terms have a strong influence on the time-mean response to an imposed forcing.

Another motivation for writing this paper was to counter some popularly held views about climate change (see, e.g., Pearce 1997; Morton 1988). First, even though recent climate change can be diagnosed in terms of known patterns of natural variability (e.g., Hurrell 1995, 1996), this does not disprove an anthropogenic cause. Second, there need be no dynamical inconsistency between anthropogenically induced climate change, and the observation of hemispheric warming trends confined to the near surface (e.g., Christy 1995).

Section 3 sets out a theoretical basis for the response of a nonlinear system with quasi-stationary regimes, to a small imposed forcing. The projection of the imposed forcing onto bases of so-called left and right singular vectors of the climate tangent dynamics is emphasized. In section 4 we discuss the evidence for the existence of quasi-stationary weather regimes, and an attempt is made to show a consistency between recent climate trends and changes in the residence frequency associated with these regimes. Supporting modeling evidence for such an association is described in section 5. The evidence that large-scale low-frequency climate may be sensitive to an imposed forcing in localized regions of space and time is given in secion 6. Section 7 gives some discussion and concluding remarks, including some pertinent suggestions for GCM evaluation.

## 2. Dynamical properties of simple nonlinear systems

### a. An analog device

In order to make more intuitive (certain aspects of) the ideas discussed in this paper, consider first a simple analog device comprising a funnel, aligned as precisely as possible above a ridge joining two channels (see Fig. 1).

A ball is thrown into the funnel and thereby drops onto the ridge joining the two channels. Under normal circumstances there is a 50% probability that the ball ends up in either one of the two cups. Some period of time after the ball has come to rest in one of the cups (greater than the time it takes to fall from the funnel into the cup), the ball is taken from the cup and the trial is repeated (many times). Due to random perturbations (arising from variations in the way the ball is dropped into the funnel, for example) we suppose that any one trial is not exactly repeatable. The resulting probability density function (PDF) of the horizontal position of the ball can be taken as bimodal with the two local density maxima corresponding to the cups. The time-mean horizontal position of the ball (the PDF mean) is directly under the funnel.

We now apply a constant weak external forcing by blowing an airstream uniformly across the whole apparatus in a given horizontal direction (e.g., with a fan). The only time that the airstream has any significant affect on the ball is when the ball is falling from the funnel to the channel ridge. As before, the ball ends up in one cup or the other, the position of the cups has not changed (so the positions of the modes of the PDF are unchanged), but there is now a higher probability that the ball will end up in one of the cups.

With this imposed forcing, and in the simplest symmetrical arangement of the channels and cups (Fig. 1 top), the ball is more likely to be found in the cup that is downstream of the funnel. In this situation, the time-mean horizontal position of the ball will move from a location at the channel ridge, to a position nearer the more probable (downstream) cup. However it is easy to deform the channels (Fig. 1 middle), so that, with the externally imposed airstream, the ball is more likely to be found in the upstream cup. In this case, the time-mean position of the ball is displaced in the opposite direction to the external forcing.

In this analog system, the response to the external forcing (the “signal”) is directly proportional to, and in the direction of, the system’s natural variability (the“noise”), irrespective of the direction of the forcing. For example, if the cups are distant 2Δ*X* apart, and the PDF is approximated by two delta functions at the cups, then the standard deviation of the unperturbed PDF system is equal to Δ*X,* and the mean response to the external forcing is Δ*p* Δ*X,* where Δ*p* is the asymmetry in the probability of the ball falling in the two cups, induced by the forcing. Hence, if the channels were lengthened so that the cups were further apart, then the time-mean position of the ball resulting from the perturbing airstream would be further deflected from the channel ridge, giving rise to an increase in signal. But the standard deviation of internal variability would also increase, giving rise to an increase in noise.

Of course, the basic system can be readily generalized with multiple regimes, as shown in the bottom panel of Fig. 1.

In this example, the externally imposed forcing is considered to be weak. If the forcing is strong, then the general character of the response may well change. With a strong enough airstream then the motion of the ball will be affected by the external forcing when it is rolling in the channels, and (e.g., in the case of Fig. 1 middle), this may be enough to inhibit the ball from reaching the upstream cup. Hence, the position of the local density maxima of the PDF may therefore start to depend on the imposed forcing once it reaches a significant amplitude. In a nonlinear system, the response will not scale linearly as the external forcing is increased to large amplitude.

### b. The Lorenz model and variants

The simple device outlined above is an analog of a rather extreme example of the inhomogeneity of the climate attractor *C* (implied by nonlinearity), with the ball representing the state vector of the dynamical system *D.* These properties can be demonstrated explicitly in simple chaotic systems such as the three-component Lorenz (1963) model whose equations of motion are

with *f*_{0} = 0. [A four-regime extension of the Lorenz model, cf. Fig. 1 bottom, is discussed in Molteni and Corti (1998).]

Figure 2 shows the impact of the additional external forcing [with *f*_{0} = 2.5 in Eq. (1)] on the PDF of the Lorenz model, projected in the *X–Y* plane for a selection of values *θ.* The dominant empirical orthogonal function (EOF) of the system points between the two regime centroids. We can define a secondary EOF to be normal to the dominant EOF in this *X–Y* plane. In the absence of any imposed forcing, the PDFs associated with the two regimes are equal. In Fig. 2 we show the effect on the PDF of forcings pointing in various directions *θ.* The first point to notice is that the phase space position of the regime centroids does not change as the forcing rotates from one direction to another. Indeed the position of the centroids is essentially unchanged if the imposed forcing is switched off altogether. This is consistent with the discussion in section 3 and the fact that the system is relatively stable in the neighborhood of the regime centroids, and relatively unstable near the origin of phase space (see Mukougawa et al. 1991).

With *θ* = 50° and 90°, the PDF associated with the upper-right regime is substantially larger than the PDF of the lower-left regime. With *θ* = 140°, the projection of the forcing in the direction of the dominant EOF is (just) pointing to the bottom left regime. However, the PDF of the upper-right regime is still larger than the PDF of the lower-left regime. In this sense, the time-mean response to the forcing in the direction of the dominant EOF, is opposite to the direction of the component of the forcing along the dominant EOF (cf. Fig. 1 middle). With *θ* = 180°, the PDF of the lower-left regime is finally larger than that of the upper-right regime.

It should be noted that the forcing’s displacement of the time-mean state away from the origin in the *X–Y* plane is not constrained to be along the diagonal joining the two regimes. To quantify this, the cosine of the angle between the diagonal and the line from the origin to the displaced time-mean state has been calculated. This can be interpreted as a correlation *C* between the time-mean response and the dominant EOF of internal variability. For the four forcing angles *θ* = 50°, 90°, 140°, and 180°, the four correlation values (to three significant figures) are *C* = 0.996, 1.00, 0.932, and 0.937. Hence, in all cases, the response is principally along the diagonal; indeed for the 90° forcing, the response is, within the accuracy of the calculation, exactly along the diagonal. As such, nonlinearity would frustrate an optimal fingerprint analysis (Hasselmann 1993) to forcings pointing in the vicinity of the 90° direction.

Now by the central limit theorem, the probability distribution of time-averaged Lorenz-model states will become Gaussian is a sufficiently long time average is applied. However, for a fixed time average, the probability distribution will be bimodal providing the state vector can be made to reside sufficiently long periods of time in one regime (yet still make chaotic regime transitions). This can be achieved by adding a variable forcing *W* to the *X* equation

(with similar equation coupled to the *Y* variable). If *X* were a stochastic white-noise variable, and *γ* = 0, then Eq. (1) would be equivalent to a Hasselmann (1976) ocean, generating a red noise (oceanic) output from a white noise (atmospheric) input. In Eq. (2), the “ocean” (*W*) is forced by a deterministic “atmosphere” (*X*), with positive feedback between *X* and *W* [cf. Eq. (1) with *W* replacing *f*_{0} cos*θ*], consistent with results from GCM integrations of the impact of midlatitude ocean variability on the atmosphere (e.g., Palmer and Sun 1985). In addition, however, the oceans may have their own states of quasi-equilibrium, for example, associated with multiple thermohaline equilibria (Stommel 1961). A simple representation of these multiple equilibria is through the nonlinear term in Eq. (2) (e.g., if *X* were stochastic, this nonlinear equation would describe random evolution in a potential well with two quasi-stable minima at *W* = ±1).

The overall effect of coupling the Lorenz model to Eq. (2) would be to to redden internal atmospheric variability, and to induce low-frequency multimodality into the *X* variability.

## 3. The response of a nonlinear system to a weak imposed forcing

In this section we review a general theoretical framework for evaluating the linearized response of a nonlinear dynamical system *D* to an imposed forcing **f**. For further details, the reader is referred to Farrell and Ioannou (1996a,b) and Palmer (1996).

The climate equations of motion, in their most general form, can be written as

where **F** is a nonlinear functional of the state vector **X**, which we assume to be (finite) *N* dimensional. Let us assume that in phase space, these equations define the“climate attractor” *C.* Now the predictability and instability properties of climate, can, for small-amplitude perturbations **x**, be derived from the linearized equation

where the Jacobian is given by

The nonlinearity in **F** implies that the Jacobian and, hence, the rate of growth of small perturbations at some point *p* ∈ *C,* varies with location on the attractor.

In addition to this inevitable consequence of nonlinearity, a second complementary aspect is emphasized in this paper. Since Eq. (3) is an *N*-dimensional set of coupled ODEs, advective nonlinearity implies the existence of quadratic terms in which the components of the state vector are multiplied together (terms involving products such as *X*_{i}*X*_{j }*i* ≠ *j*). The existence of such terms means that, in general, **J** will not be symmetric (i.e., *J*_{ij} ≠ *J*_{ji} for at least some *ij*). These two aspects of nonlinearity can be illustrated in the simple example

where

is clearly dependent on *X*_{1} and *X*_{2}, and for *X*_{1} ≠ *X*_{2}, *J*_{12} ≠ *J*_{21}. More generally, for this type of nonlinearity, **J** will be non-normal, that is, **J*****J** ≠ **JJ***, where **J*** denotes the adjoint of **J**. (Non-normality rules out the possibility of antisymmetric or complex-Hermitian operators.) The consequences of non-normality are discussed below.

Now Eq. (4) can be written in the integral form

where **L**(*t*_{1}, *t*_{0}) is the forward tangent propagator. This propagator can be related to the Jacobian by splitting the nonlinear trajectory **X** into many short quasi-stationary segments and writing, for each segment *t*_{i} < *t* < *t*_{j},

We now perturb the nonlinear system in Eq. (3) with some weak imposed forcing **f**(*t*). In this case, the linearized dynamical Eq. (4) becomes

Using the tangent propagator, the solution to Eq. (9) can be written over the finite interval [*t*_{1}, *t*_{0}] as

If **f** is time independent, then

where

The system’s response to **f** depends on the projection of **f** onto local phase-space instabilities on *C.* To see this more explicitly, let us perform a standard singular vector decomposition (e.g., Strang 1986) on ℒ, so that

where **U** is the matrix whose columns contain the left singular vectors **u**_{i} of ℒ, and **V** is the matrix whose columns contain the right singular vectors **v**_{i} of ℒ, all relative to (and dependent on) the time interval [*t*_{1}, *t*_{0}]. The right singular vectors **v**_{i} satisfy the eigenvector equation

The left singular vectors **u**_{i} satisfy the corresponding eigenvector equation

Finally, in Eq. (13), **Σ** is a diagonal matrix whose real elements, the singular values Σ_{ii} = *σ*_{i} > 0, are ordered so that *σ*_{1} ≥ *σ*_{2} ≥ . . . ≥ *σ*_{N}.

Now if **J** is non-normal, so is ℒ. Hence, by definition, the right and left singular vectors are not equal to one another. On the other hand, since the compound operator ℒ*ℒ is necessarily normal (irrespective of whether ℒ itself is normal) then the left and right singular vectors each form an orthogonal set. We assume that this set is complete in the *N*-dimensional phase space, so that the left and right singular vectors can each be chosen to form an orthonormal basis. The left singular vectors singular vectors **u**_{i} will approximate the Lyapunov vectors of the dynamical system *D* (see Legras and Vautard 1995; Smith 1996; Toth and Kalnay 1993 in addition to the references at the beginning of this section). In physical terms, the dominant Lyapunov vectors describe the principal instabilities of the system.

The relationship between the dominant Lyapunov vectors and the regime centroids is not immediate. Broadly speaking, one can only expect a correspondence providing the fastest Lyapunov timescale is equivalent to a typical regime transition timescale. Since the latter is on the order of days, a correspondence will only be possible for atmospheric models with suitable balance constraints (e.g., quasigeostrophic balance) to filter very rapid subsynoptic-scale (e.g., convective) instabilities.

It can be noted that under the conditions of the Fluctuation–Dissipation theorem (Leith 1973), the response to an imposed forcing is determined by the lagged covariance matrix of the undisturbed system. This theorem is consistent with (10)–(13) when the singular vector pairs that have large singular values correspond to components of the lagged covariance matrix with long lags. The precise conditions under which the Fluctuation–Dissipation theorem holds (which require, among other things, the climate attractor to have a Gaussian PDF) are not met in the global climate system. Nevertheless, the general perspective put forward here is qualititatively consistent with the statement of the Fluctuation–Dissipation theorem, even though departures from strict Gaussianity are emphasized in this paper.

Now let 〈 · , · 〉 denote the inner product with respect to which the adjoint operator is defined, and assume that 〈**f**, **f**〉 = 1. Let us represent the nonlinear quasi-stationary regime anomaly fields by the vectors **E**_{i}. Then the normalized forcing, which will optimally excite **E**_{i}, is one in which 〈ℒ**f**, **E**〉 is maximized. This optimal forcing, or sensitivity pattern is therefore defined by

If we apply the singular vector decomposition [Eq. (13)], then

Now let us assume that our model is sufficiently meteorologically balanced that the projection of **E** onto the singular vectors **U** is either weighted toward the leading singular vectors (as is the case in the quasigeostrophic model used in this paper, see section 6 below), or is, at worst, uniform among the spectrum of singular vectors. Then, from Eq. (17), the sensitivity vector will be weighted, through the matrix **Σ**, toward the dominant right singular vectors. As will be shown explicitly, the sensitivity vector can be rather localized in space. This is because, as shown by Buizza and Palmer (1995), the evolution of the dominant singular vectors (from initial to optimization time) involves an upscale transfer of energy. Sensitivity vectors have proved an extremely useful tool in the diagnosis of forecast error for weather prediction (Rabier et al. 1995).

In view of the weighting of the sensitivity vector toward the dominant right singular vectors, its practical importance depends on the projection of the imposed forcing **f** onto the right singular vectors. We will show below that the projection of a planetary-scale forcing **f** in a quasigeostrophic model projects reasonably uniformly onto a basis of right singular vectors. For such a projection, the sensitivity vector defines the most important component of **f** as far as the excitation of the regime centroids is concerned.

Because of the nonlinearity in **F**, the dominant singular values will vary in magnitude around the chaotic attractor. Regions where these values are relatively large in magnitude denote regions of local instability and will therefore not be associated with quasi-stationary local density maxima. On the other hand, regions where the dominant singlar values are relatively small in magnitude need not necessarily be associated with regimes; for example, there may be regions of phase space through which the state vector is in stable transition between regimes. A study of variations in singular values and quasi-stationary states has been described by Yamane and Yoden (1997).

If the state vector is relatively insensitive to **f** in the vicinity of a regime centroid, then the (geographical) structure of the regime should be relatively unaffected by the imposed forcing. By contrast, we can expect **f** to have a much stronger influence on the PDFs associated with these regimes since these PDFs can be affected by parts of attractor remote from the regimes themselves.

We have emphasized two complementary aspects of nonlinearity. First, that **J** is a function of position on the attractor, and second, that **J** is non-normal. These considerations lead us to propose the following nonlinear perspective based on a system whose variability is determined principally by two or more quasi-stationary regions between which the system makes frequent transitions.

The response of climate to a small imposed forcing

*f*will be manifest primarily in terms of changes to the PDFs associated with the quasi-stationary regimes.The geographical structures of these regimes will be relatively insensitive to

*f.*The influence of

*f*may be strongest in rather localized regions of space and time (corresponding to the regime sensitivity patterns).

Evidence for this perspective will be discussed in sections below.

## 4. Quasi-stationary regimes and climate trends

### a. Evidence for regimes

The existence of quasi-stationary regimes have been found explicitly both in intermediate models of the atmospheric global circulation (e.g., Reinhold and Pierrehumbert 1982; Legras and Ghil 1985), and in GCMs (Haines and Hannachi 1995; Hannachi 1997; Robertson et al. 1998). (In the context of the dynamical perspective developed here, the realism of a GCM’s representation of regime activity is an important issue in the validation of the model, see section 7.)

Can such quasi-stationary regimes be diagnosed from climate data? From a synoptic point of view, the existence of persistent anomalies (with transition times much shorter than residence times) has been well established (e.g., Oerlemans 1978; Dole and Gordon 1983;Yang and Reinhold 1991). Corresponding to this, there seems little doubt that weather regimes can be defined on a sectorial basis, using either cluster analysis or some method that tests for quasi-stationarity or non-normality of a corresponding distribution function (Vautard 1990;Kimoto and Ghil 1993b; Cheng and Wallace 1993; Robertson et al. 1998; Smyth et al. 1998, manuscript submitted to *J. Atmos. Sci.*). The corresponding regime centroids bear resemblence to teleconnection patterns or persistent anomaly fields (Wallace and Gutzler 1981; Dole and Gordon 1983). For example, anomaly fields for weather regimes defined over the Pacific sector are highly correlated with the Pacific–North American (PNA) pattern. It can be noted, on the other hand, that the high and low phase of the North Atlantic oscillation (NAO) pattern is not so clearly defined in terms of individual Atlantic sector regimes.

The existence of hemispheric-scale regimes is more controversial. For example, while early studies of bimodality on hemispheric-scale planetary-wave amplitude distributions were claimed by Hansen and Sutera (1986), the statistical significance of these results are not universally accepted (Nitsche et al. 1994). Moreover, doubts have been raised as to the physical coherence of hemispheric regimes. For example, Wallace (1996) notes that the PNA and NAO patterns can fluctuate independently for much of the time.

On the other hand, cluster analysis studies (Mo and Ghil 1988; Molteni et al. 1990; Kimoto and Ghil 1993a) focusing on a phase space spanned by hemispheric EOFs have produced reasonably robust hemispheric-scale regimes (see Fig. 3).

Moreover, there is some evidence for hemispheric coherence. For example Gutzler et al. (1988) show non-negligible correlations between lower-tropospheric temperature over the whole hemisphere, and a regionally defined PNA index. In addition, there is evidence of an (interannual-timescale) correlation between satellite-derived snow cover over North America and Eurasia (Gutzler and Rosen 1992; Walland and Simmonds 1997).

It is possible that this conflicting evidence for the existence of hemispheric regimes can be reconciled if Pacific and Atlantic sector variability, although predominantly quasi-independent, can become partially synchronized from time to time. The notion of atmospheric chaotic synchronization has been discussed in the context of interhemispheric synchronization (Duane 1997), and may therefore also apply to intersector synchronization. It should be noted that a linear correlation analysis may not be the most appropriate means of assessing the existence of possible partial synchronization.

It can also be noted (Molteni and Corti 1998) that unambiguous evidence for regime structure within a reduced-dimensional phase space may be difficult if the analysis includes periods when the state vector is mainly fluctuating in a subspace orthogonal to the reduced space.

The basic nonlinear paradigm put forward in this paper may be more appropriate to the PDFs and centroids of the sectorial regimes, if, as appears likely, these are more fundamental to the dynamics of the atmosphere than the hemispheric regimes. On the other hand, the discussion in this paper is in the context of observations of trends in hemispheric-mean surface temperature, and as such we attempt to link such trends to changes in the residence frequency of the hemispheric regimes shown in Fig. 3. It is hoped in a future paper to apply a more quantitative application of this methodology on the basis of a sectorial regime analysis.

It should be noted that although only two clusters are shown in Fig. 3, in practice, both sectorial and hemispheric regime analyses define multiple centroids; for quantitative purposes, it is impossible to reduce the number of clusters to just two, as in the examples in section 2.

The extent to which interannual and interdecadal fluctuations in the atmospheric flow has multimodal behavior, may depend on the nature of the coupling to the oceans. While atmospheric coupling to a Hasselmann ocean [cf. Eq. (2) with *γ* = 0] may be sufficient to allow multimodal behavior on interannual timescales, the existence of multiple equilibrium steady states associated with the thermohaline circulation (Stommel 1961) provides a possible mechanism to generate significant atmospheric multimodality on decadal and longer timescales [Eq. (2) with *γ* ≠ 0]. The clear existence of multimodality on these timescales, from either observational data, or from GCM data, remains to be established.

Even with such ocean coupling, the spatial structure of interannual and interdecadal regimes can still be expected to be similar to that determined from the intraseasonal analyses referenced above (the SST derived from the coupling can be thought of as a time-varying weak forcing on the atmosphere). In this respect, it can be noted that the dominant EOFs of seasonal (Wallace 1996) and decadal (Trenberth and Hurrell 1994) average winter states correspond well with the regime structures shown in Fig. 3.

In the Tropics, the dominant pattern of interannual variability is associated with El Niño. During El Niño years, the Tropics as a whole is warm, and, moreover, the El Niño has a tendency to increase the frequency of cluster 2 events (see below).

Although on timescales of a few months, much of the El Niño can be understood in terms of linear dynamics, there is evidence of nonlinear chaotic dynamical structure (e.g., Münnich et al. 1991). A nonlinear dynamical systems analysis of the El Niño phenomenon (Wang et al. 1999) has suggested that in certain circumstances, the El Niño system may possess two quasi-equilibrium (warm and cold) states, with transitions between the states triggered by stochastic forcing [cf. Eq. (2) with *γ* ≠ 0]. The possible bimodal nature of El Niño was first noted by Wyrtki (1982), though again this is controversial. For example, although Penland and Sardeshmukh (1995) did not find evidence for regime dynamics, recent analysis of the distribution of SST anomalies in the so-called Niño 3.4 region (Trenberth 1997) gives observational support for possible bimodality.

However, irrespective of the evidence for multimodality, the importance of nonlinearity in the tropical Pacific has been shown explicitly in a GCM context by Stockdale et al. (1993). An ocean GCM was run over 9 yr forced first with monthly mean wind stresses, and second with the 9-yr seasonal mean of the wind stresses. The tropical west Pacific was 0.7°C cooler in the first integration than in the second; the differences arising from nonlinear rectification. Such an SST difference field could generate a substantial impact on the PDF of the northern hemisphere winter regimes, leading to signficant decadal average hemispheric-mean surface temperature differences.

### b. Regime frequency and climate trends

Consistent with Gutzler et al.’s (1988) analysis, a regression between hemispheric mean temperature and 500-hPa geopotential height gives rise to a pattern that is well correlated with the cluster 2 pattern in Fig. 3 (Wallace et al. 1996; see also Hurrell 1996). The corresponding surface temperature regression pattern shows, to first approximation, cooling over the oceans, and warming over land; the so-called Cold Ocean–Warm Land (COWL) pattern. Hence, the cluster 2 regime centroid can be thought of as a midtropospheric representation of the COWL pattern, with an anomalously warm hemispheric-mean surface temperature.

Now Fig. 4 shows the geographical distribution of the change in surface temperature between the 1950s and the 1980s, as reported in IPCC (1992). The increase in surface temperature is largest in winter, where it is concentrated over the North America and Northern Eurasia, and is partially offset by cooling over the North Atlantic and the North Pacific. The change in surface temperature in Fig. 4 is therefore consistent with an increase in the PDF associated with cluster 2 in Fig. 3, (and a decrease in the PDF associated with cluster 5).

If the horizontal structure of recent climate trends is consistent with the horizontal structure of a dominant quasi-stationary regimes, then the vertical structure of the same climate trend should similarly be correlated with the vertical structure associated with this regime.

As mentioned in the introduction, popular articles have drawn attention to the IPCC (1996) analysis that while there are significant trends in global temperature based on surface measurements, the corresponding trends from midtropospheric radiance measurements are much weaker. [Pearce (1997) cites these contrasting sets of observations as an apparent focus for “Greenhouse Wars” between climate scientists, as to whether anthropogenic global warming is really occurring.]

As Hurrell and Trenberth (1996) have noted, satellite and surface measurements of global temperature change give different perspectives on the same events. In the context of the dynamical perspective put forward in this paper, the contrasting trends from the two sets of observations are entirely consistent with anthropogenic forcing. The hemispheric warmth of the cluster 2 regime at the surface is a consequence of its spatial structure, with positive height centers over land areas, negative height centers over oceanic areas. As noted above (and in Hurrell and Trenberth 1996), variability in surface temperature is large over land and small over the ocean, due, in turn, to the relatively small heat capacity of the active land, and the ability of the ocean to mix down surface anomalies.

According to the perspective developed here, the observed climate change in the Northern Hemisphere can be interpreted as an increase in the PDF of the climate attractor associated with the cluster 2 regime. Hence, the vertical structure of the observed hemispheric warming should be correlated with the vertical structure of the weather regimes, and less the vertical structure of anthropogenic forcing (e.g., toward a one-dimensional radiative equilibrium temperature profile associated with an increased greenhouse effect). The hemispheric-mean temperature anomaly associated with the cluster 2 regime is primarily a surface phenomenon. In the free troposphere, the temperature anomalies are determined by the equivalent barotropic nature of the regime structure, and are similar in size over the ocean and over the land. As such, the cluster 2 pattern will not be associated with substantial hemispheric-mean temperature anomalies in the free troposphere.

Hence, if the observed increase in hemispheric warming is an anthropogenically induced increase in the PDF of the cluster 2 regime, then such an anthropogenic influence should have relatively little impact on hemispheric-mean temperatures in the free troposphere.

In the IPCC (1990) report, it was noted that in the upper troposphere (300–100 hPa), observations reported by Angell (1988) showed a rather steady decline in temperature, in general disagreement with many model simulations that warm at these levels when the concentration of greenhouse gases is increased. (The simulations show mean cooling, only in the stratosphere.) Again, these observations can be readily understood using the dynamical perspective developed here. The analysis of Angell was based on radiosonde observations, and therefore describes temperature trends mainly over land. Now, the geopotential height anomalies of the cluster 2 regime are largely positive over land, and the vertical structure of such large-scale anomalies are mainly equivalent barotropic with maximum geopotential amplitude at the equivalent barotropic level of about 300 hPa. Hence, over land, the equivalent barotropic structure will ensure positive temperature anomalies below about 300 hPa and negative temperature anomalies above about 300 hPa.

Of course it can be asked why many GCMs forced with enhanced CO_{2} do not replicate the observed vertical structure. Of course it is very likely that effects from aerosol and ozone forcing cannot be neglected. Also vertical resolution near the tropopause may be inadequate in many models. However, it is also possible that (due to model systematic error associated with zonalization), both regime structure and corresponding regime sensitivity, may be quantitatively weaker in some of the current generation of models, than in the real atmosphere.

## 5. Evidence for the nonlinear perspective from GCM integrations

In this section, evidence from atmospheric GCM integrations with prescribed forcing of a response consistent with a change in regime frequency is discussed. In most (but not all) of these cases, the prescribed forcing derives from the anomalies in the lower boundary.

The ability of tropical Pacific SST anomalies to influence the PDFs of the northern hemisphere quasi-stationary winter regimes on interannual timescales was first discussed in Molteni et al. (1993). Figure 5 illustrates an analysis of two (three-member) 120-day ensembles of the ECMWF atmosphere model. The first ensemble was run with observed SSTs for 1986/87 (when tropical Pacific SST was anomalously warm), the second ensemble with observed SSTs for 1988/89 (when tropical Pacific SST was anomalously cold). The first panel shows that the PDF associated with cluster 2 (see Fig. 3) is larger in 1986/87 compared with 1988/89, while the PDF asssociated with cluster 5 is larger in 1988/89 than in 1986/87. The time-mean ensemble-mean difference patterns shown in the lower panels in Fig. 5 are similar to the difference between the cluster 2 and 5 patterns themselves, consistent with the dynamical perspective of this paper.

On the other hand, it is simplistic to imagine that all such interannual variability can be represented in terms the frequency change solely of these two hemispheric regimes. However, the impact of interannual SST fluctuations on sectorial regime frequencies has been discussed elsewhere, both over Europe (Fraedrich 1990) and the United States (Robertson and Ghil 1999). Such sectorial analyses lead to practical probability forecasting products.

The ability of tropical Pacific SST anomalies alone to generate the hemispheric-scale COWL pattern on decadal timescales has been conclusively demonstrated by Kumar and Hoerling (1998). Specifically, these authors have simulated the observed decadal warming of Northern Hemisphere land surface temperatures in the years from 1978 to 1994, using observed SST anomalies in the tropical Pacific Ocean only. Clement et al. (1998) have taken this a stage further. Observed SSTs in the Pacific only have been used to force an atmospheric GCM coupled in middle latitudes to a mixed layer ocean model. The resulting midlatitude SST anomalies showed broad cooling in the North Atlantic and North Pacific, consistent with the COWL pattern over the northern oceans.

A sectorial analysis of changes in regime frequency on decadal timescales over the Atlantic, as a result of imposed SST anomalies has been recently carried out by Robertson et al (1998). Specifically, weather regimes over the North Atlantic sector were first determined from a 100-yr control integration of the Max Planck Institute ECHAM3 model, with seasonally varying climatological SST. The resulting six regimes were found to be quite realistic. A second experiment with an observed decadal SST anomaly over the North Atlantic was then analyzed in a similar way. The SST anomaly was found to change the frequency of certain weather regimes by about one standard deviation. On the other hand the regime patterns in the anomaly experiment were not statistically significantly different from the regime patterns in the control experiment. Robertson et al (1998) construct a response pattern to the SST anomaly by computing the mean of the six weather regimes, weighted according to the changes in the residence frequency for the regimes. The authors note that this gives a more robust picture of the response to the imposed SST anomaly than that determined by the more conventional time-mean response.

As far as the author is aware, an analysis of integrations with and without doubled CO_{2} similar to that performed by Robertson et al (1998) has not yet been performed. Perhaps the closest is that by Timmermann et al. (1998) who describe results from an integration of the ECHAM coupled model forced with monotonically increasing atmospheric CO_{2} concentrations. The EOF of surface temperature that best described the response to an increase in CO_{2} was found to be strongly spatially correlated with a certain decadal mode of internal variability of the coupled ocean–atmosphere system. Although much of the air–sea coupling in this model occurs over the North Atlantic involving the thermohaline circulation, the atmospheric component of the decadal mode was described reasonably well by the hemispheric-scale regimes shown in Fig. 3. In this sense the results of Timmerman et al. (1998) support the general notion put forward here, that the signal arising from anthropogenic climate change has a substantial component associated with a change in the PDFs of the dominant natural modes of variability.

In section 4 it was suggested that nonlinearity associated with air–sea coupling in the tropical Pacific, and the associated impact of El Niño on global temperature indicates that it is of importance to be able to estimate the impact of an externally imposed forcing on the PDF associated with El Niño. Consistent with this, E. Schneider (1997 personal communication) found that the impact of doubling CO_{2} on global-mean temperature in an atmosphere–mixed layer-ocean model reduced substantially when the SSTs in the eastern equatorial Pacific cold tongue were held constant at the control simulation values. (The reduction was out of proportion to the size of the area in which SST was held constant.)

## 6. Sensitivity of regime PDFs to an imposed forcing

In this section we study the third aspect of the nonlinear perspective as outlined at the end of section 3 that the climate may be particularly sensitive to an imposed forcing, in localized regions of space and time.

### a. Sensitivity in a baroclinic time-varying flow

Corti and Palmer (1997) have calculated the sensitivity of the PNA pattern to initial perturbations 5 days earlier, in a T21 three-level quasigeostrophic model (Marshall and Molteni 1993). Similar calculations for atmospheric blocking patterns have been reported by Oortwijn and Barkmeijer (1995). The sensitivity pattern is defined as (see section 3)

where **E**_{PNA} is a fixed pattern, matching the observed PNA pattern, but determined from an EOF analysis of a long integration of the quasigeostrophic model. In the quasigeostrophic model, the large-scale regimes project well onto the dominant left singular vectors **u**_{i} [this has been shown explicity by Oortwijn (1998) for the Euro-Atlantic blocking regime]. As such **S**(*t*) will have maximum projection onto the corresponding leading right singular vectors **v**_{i}.

Figure 6 shows an example of the time variation of ‖**S**(*t*)‖ together with the time variation of the model’s PNA index (the projection of the state vector onto **E**_{PNA}). As discussed above, the time variation of **S**(*t*) is consistent with the nonlinearity in **F**. Note also that ‖**S**(*t*)‖ varies on rather shorter timescales than does **E**_{PNA}. [It can be seen that ‖**S**(*t*)‖ is not well correlated with the PNA index. As can be seen periods of large ‖**S**(*t*)‖ are rather localized in time. In fact, as discussed in Corti and Palmer (1997), ‖**S**(*t*)‖ is correlated with a secondary EOF orthogonal to **E**_{PNA}. This is consistent with what is found in the Lorenz model; variations in the sensitivity pattern associated with the first Lorenz-model EOF are correlated with the principal components of the second EOF that points in the *Z* direction, and along which the singular values have their greatest variation.]

An example of the PNA sensitivity vector is shown in Fig. 7. Like the adjoint normal mode from the barotropic model (see Palmer 1993), the sensitivity vector is relatively localized between 90°E and the date line. Unlike the barotropic adjoint, the sensitivity vector has a baroclinic tilt (not shown), with most amplitude in the midlower troposphere. The reason for this can be understood in terms of wave-action dynamics (see Buizza and Palmer 1995).

The relative efficiency with which **E**_{PNA} can be forced using an averaged sensitivity pattern is demonstrated in Fig. 8. The top panel of this figure shows the PDF of the quasigeostrophic model PNA index based on a control integration (solid line). Also shown is the PDF of the same index when the model is forced with ± a weak imposed time invariant forcing proportional to **E**_{PNA}. By contrast the bottom panel of Fig. 8 shows the control and perturbed PDFs, where the forcing has the same amplitude as above, but is proportional to a 2000 case-average sensitivity field. It can be seen that the perturbed PDFs are considerably further displaced from the control PDF when the forcing is proportional to the case-average sensitivity pattern, than when it is proportional to **E**_{PNA}.

It can be noted in passing that the PDF of the PNA pattern is noticeably broader when the imposed forcing is toward the direction of negative PNA index. This is consistent with the fact that the negative PNA regimes have higher internal barotropically induced variability, than the positive PNA regimes (Palmer 1988; Molteni and Palmer 1993).

Of course, the choice of a 5-day period for these sensitivity calculations is ad hoc. However, the longer the period over which the sensitivity is computed, then the more the sensitivity vector will become dependent on flow details, and hence the more the spatial structure of the sensitivity vector will vary from case to case. As such, the case-averaged sensitivity vector (giving the time-mean sensitivity forcing) will become less and less coherent as the time over which the sensitivity is computed gets longer and longer. On the other hand, over very short periods, the sensitivity, while very coherent from case to case, will become small in amplitude. Overall, this suggests that an optimum period may exist; that is, it may be possible to estimate an optimal sensitivity to forcing (with respect to time-dependent basic state dynamics) without explicit reference to an optimization time.

At this stage, it is not possible to be able to claim any quantitative significance to the sensitivity calculations (at least for climate prediction). On the other hand, the calculations illustrate (in the context of a reasonably realistic atmospheric model) the third component of the nonlinear perspective: that the climate (attractor) is sensitive to an imposed forcing, in localized regions of space and time.

### b. Projection of a planetary-scale forcing onto a basis of singular vectors

The entire spectrum of singular vectors of the T21L3 quasigeostrophic model (used above) has been computed by Reynolds and Palmer (1998). Figure 9 shows (48-h) singular values, *σ*_{i} from bases of (1449) singular vectors of **L** (averaged over six different cases taken from ECMWF analyses during the northern winter season). It can be seen that there is a substantial slope to this singular-value curve, with *σ*_{1} being over three times larger than *σ*_{100}, and nine times larger than (the neutral) *σ*_{500}.

Examples of a leading, neutral, and decaying right singular vector from the quasigeostrophic model are shown in Reynolds and Palmer (1998). The dominant right singular vector is similar to the PNA sensitivity perturbation shown in Fig. 7, that is, showing strongly localized structure over the Asian landmass. By contrast, the neutral singular vector is less localized with a maximum in the upper troposphere. The trailing singular vector is also localized over the Asian continent and the west Pacific, though with scales larger than those of the leading singular vector.

We study the projection of the gravest possible spectral component of the quasigeostrophic model state vector (considered to represent a possible planetary-scale forcing) onto the entire basis of quasigeostrophic right singular vectors. Specifically we choose for **f** the spectral component of quasigeostrophic streamfunction that has wavenumber zero in the zonal direction (that is, zonally symmetric) and wavenumber 1 in the meridional direction. The top panel of Fig. 10 shows the corresponding projection coefficients *α*_{i}, where the *i*th data point represents the value of *α*_{i} averaged over the six chosen cases. It can be seen that the *α*_{i} do not vary very smoothly from one singular vector to the next. Nevertheless some nonmonotonic low-frequency variation emerges, with relative maxima in the unstable and neutral subspaces.

However, this low-frequency variation is rather small compared with the variation in the *σ*_{i}. In particular, the bottom panel of Fig. 10 shows the variation of the product *σ*_{i}*α*_{i}. The monotonic slope of the singular value curve overwhelms the variation between the *α*_{i}, so that the “low-frequency” component of the slope of the product *σ*_{i}*α*_{i} curve is monotonic with largest values in the unstable subspace. In this sense, the the planetary-scale forcing **f** projects reasonably uniformly onto the basis of right singular vectors. As given by the analysis in section 3, the response to such a forcing will be mainly determined by the dominant left singular vectors, and hence the regime anomalies.

(Note that in this analysis, there is nothing special about the choice of a 48-h optimization time. If the prescribed forcing pattern **f** projects quasi-uniformly onto one arbitrarily chosen singular vector basis, it is likely to project quasi-uniformly onto all chosen singular vector bases.)

### c. Sensitivity of interanual and interdecadal modes of variability

As with weather regimes, the El Niño PDF can be most readily changed by an imposed forcing that projects onto the El Niño sensitivity pattern. Right singular vectors of El Niño have been computed in intermediate coupled models, for example, in Chen et al. (1997). Typically, the SST patterns associated with these singular vectors have a dipole structure, for example, with cold SSTs in the western Pacific and warm SSTs in the eastern Pacific. The corresponding left singular vectors have the structure of a mature warm El Niño event. The associated singular values show marked variability with time of year, and the state of El Niño in the basic-state trajectory. More pertinent in the present context, Moore and Kleeman (1999) have computed the stochastic optimals (dominant singular vectors of ℒ) for the El Niño event. The SST structures are not disimilar to Chen et al.’s right singular vectors.

On the decadal timescale, some work on the optimal forcing of the thermohaline circulation has been reported by Lohmann and Schneider (1998), who have made a singular vector analysis of the nonlinear box model of Stommel (1961). Lohmann and Schneider (1998) find that the leading right singular vector has substantial projection onto high-latitude haline fields.

It would clearly be more satisfactory if it were possible to derive some generic sensitivity fields in which climate fluctuations on all timescales could be treated in a unified consistent way (rather than the piecemeal study of day–day, year–year and decade–decade sensitivity fields as described here). For example, in computing El Niño singular vectors, much of the relatively fast synoptic timescale atmospheric processes have been filtered from the equations of motion (just as subsynoptic timescales are filtered in the quasigeostrophic model). Similarly, the model used for computing thermohaline singular vectors are filtered for interannual fluctuations. Unfortunately, it is not clear how to estimate El Niño or thermohaline sensitivity including the faster timescales in the equations of motion, not least as there is no unambiguous way to treat the nonlinear saturation of these fast timescale processes over long periods.

For example, consider three imposed time-invariant normalized forcings, the first given by an average 5-day sensitivity vector of the NAO, the second given by a thermohaline sensitivity vector, the third given by some random stochastic noise. Which of the forcings would have the largest impact on a very long time-average NAO index? While it would appear unlikely to be the third, it is difficult a priori to choose between the first and the second. However, in principal this sort of question can be studied by applying the three forcings (the former two estimated from linear models) into the equations of motion of a (nonlinear) GCM. If there was a clear-cut answer, it would help clarify the role that this type of linearized analysis can play in studying climate sensitivity.

Nevertheless, irrespective of these considerations, the sensitivity estimates (on all the timescales considered) are consistent with the third component of the nonlinear perspective outlined in section 6 above, that is, large-scale climate anomalies can be sensitive to forcing in relatively localized regions of space and time.

From a modeling point of view, the broad conclusion to be drawn from this conclusion is that it may not be permissable to ignore errors in simulating climate that occur on space and time scales short compared with the scales of interest.

## 7. Discussion

In this paper, we have considered the response of climate, considered as a nonlinear dynamical system *D,* to a weak imposed forcing **f**, representing the effects of anthropogenic increases in CO_{2}.

Of course, it is not immediate how to translate an increase in CO_{2} to some specified **f**. For example, as far as a GCM is concerned, the CO_{2} mixing ratio is a model parameter whose value is specifiable at the start of an integration. Obviously, the most direct impact of an increase in CO_{2} is through the radiative (greenhouse) forcing of the atmosphere. As such, **f** could be taken to denote a global-scale anomalous forcing in the thermodynamic equation (toward warmer temperature). Such a (purely thermodynamic) forcing would have no direct projection onto (for example) the streamfunction fields associated with the illustrated singular vector and sensitivity patterns. On the other hand, on very short timescales a pure thermodynamic forcing will induce perturbations in divergence and vorticity as the system adjusts towards a balanced state. As such the impact of enhanced CO_{2} can be thought of in terms of some imposed balanced forcing **f**_{T,D,ζ}, which would then have some direct projection onto the sensitivity patterns discussed above.

Now, in developing the dynamical perspective, it is not so important to have a clear picture on the precise form of **f**_{T,D,ζ}. Rather the perspective is based on gross inhomogeneities in the climate attractor, and on general characteristics of the projection of **f**_{T,D,ζ} onto local instabilities of the attractor. In summary, the basic constituents of the dynamical perspective that we are trying to apply to the climate change problem can be summarized as follows.

We consider a system

*D*whose variability is determined principally by transitions between quasi-stationary states.We consider a weak imposed forcing

**f**on*D*that projects reasonably uniformly onto any set right singular vectors of*D*’s tangent propagator.

The basic consequences of this can be outlined as follows.

The system responds to

**f**primarily through a change in the PDF of the principal components associated with*D*’s quasi-stationary states.The physical structure of the dominant EOFs is relatively unaffected by the presence of the imposed forcing.

Even though the quasi-stationary states may be large scale, changes in the PDFs associated with these states may be particularly sensitive to

**f**in localized regions of space and time.

While the focus in this paper has been principally on hemispheric-mean climate change, the perspective is probably of more fundamental relevance to regional climate change. This is because most of the current evidence points toward the nonlinear regimes being fundamentally sectorial in nature, though with partial synchronization. One can note that the perspective may also be relevant on a regional basis in the Tropics; for example if Indian monsoon rainfall is broadly determined by the PDFs of the active and break monsoon regimes (associated with two quasi-stationary ITCZ positions; Webster et al. 1998), then the impact of enhanced CO_{2} on the monsoon climate depends on CO_{2}-induced changes to the PDF of these regimes.

None of the analysis above suggests that climate prediction using GCMs is impossible. However, the analysis does highlight the potential importance of spatially (and temporally) localized model errors. Errors in the representation of physical processes in areas (and times) where the regime sensitivity patterns have large amplitude, could have significant impact on predictions of hemispheric mean temperature. It should be noted that regional errors in the simulation of heat fluxes at the air–sea interface by coupled GCMs can exceed the 4 W m^{−2} associated with the direct effect of doubled CO_{2}, by over an order of magnitude.

The analysis above suggests a number of tests for the evaluation of climate models.

The regime structure of a GCM should be similar to that of the real atmosphere. According to the perspective in this paper, incorrect regime structure would imply an erroneous time-mean response to enhanced CO

_{2}. It is recommended that sectorial regime diagnosis techniques become standardized, in order that climate model regime intercomparisons can take place.For those models that have realistic sectorial regime structure, it is proposed that tests be performed to estimate whether the imposed increase in CO

_{2}has had a larger impact on the regime PDFs than on the regime structure. The robustness of the response to increased CO_{2}as determined by the changes in regime frequency should be compared with the full time-mean response. In this respect, determining the climate change response from only those regimes in which regime frequency has changed significantly, could provide a “fingerprint” method consistent with the nonlinear perspective in this paper.Adjoint sensitivity analysis should be performed on climate-model fields. It is possible that the zonalization and loss of planetary-wave amplitude, characteristic of many GCM integrations, would lead to less sensitivity compared with the real atmosphere. In such a situation it is possible that the GCM response to anthropogenic forcing would have a weaker projection onto the model’s patterns of internal variability than in the real world.

It is essential for a GCM to be able to simulate the spectrum of SST variability in the tropical Pacific; not least the tropical Pacific appears as a sensitive area for forcing changes in COWL-pattern frequency. In order to simulate such coupled ocean–atmosphere dynamical processes accurately, the meridional resolution of such a GCM must be a fraction of a degree. As noted in IPCC (1996), no coupled model used for projecting the response to enhanced CO

_{2}had (at the time of writing) such resolution. From the perspective of this paper, this must be viewed as a significant shortcoming of published climate change simulations.

It should be noted that adequate ocean resolution is not a sufficient requirement for modeling accurately the dynamics associated with the Pacific sensitivity region. For example, nonlinear moist processes in the atmosphere are also known to play a central role in determining SST in the tropical Pacific. One pertinent means of evaluating a climate model is to test its skill in forecasting interannual variability in the tropical Pacific region. As noted in IPCC (1990): “confidence in a model used for climate simulation will be increased if the same model is succesful when used in a forecasting mode.” The worldwide upsurge in activity in the field of seasonal forecasting bodes well for this type of model evaluation.

## Acknowledgments

## REFERENCES

## Footnotes

*Corresponding author address:* Dr. Timothy N. Palmer, ECMWF, Shinfield Park, Reading RG2 9AX, United Kingdom.

Email: t.palmer@ecmwf.int