## Abstract

A 62-member ensemble of coupled general circulation model (GCM) simulations of the years 1940–2080, including the effects of projected greenhouse gas increases, is examined. The focus is on the interplay between the trend in the Northern Hemisphere December–February (DJF) mean state and the intrinsic modes of variability of the model atmosphere as given by the upper-tropospheric meridional wind. The structure of the leading modes and the trend are similar. Two commonly proposed explanations for this similarity are considered.

Several results suggest that this similarity in most respects is consistent with an explanation involving *patterns* that result from the model dynamics being well approximated by a linear system. Specifically, the leading intrinsic modes are similar to the leading modes of a stochastic model linearized about the mean state of the GCM atmosphere, trends in GCM tropical precipitation appear to excite the leading linear pattern, and the probability density functions (PDFs) of prominent circulation patterns are quasi-Gaussian.

There are, on the other hand, some subtle indications that an explanation for the similarity involving preferred *states* (which necessarily result from nonlinear influences) has some relevance. For example, though unimodal, PDFs of prominent patterns have departures from Gaussianity that are suggestive of a mixture of two Gaussian components. And there is some evidence of a shift in probability between the two components as the climate changes. Interestingly, contrary to the most prominent theory of the influence of nonlinearly produced preferred states on climate change, the centroids of the components also change as the climate changes. This modification of the system’s preferred states corresponds to a change in the structure of its dominant patterns. The change in pattern structure is reproduced by the linear stochastic model when its basic state is modified to correspond to the trend in the general circulation model’s mean atmospheric state. Thus, there is a two-way interaction between the trend and the modes of variability.

## 1. Introduction

A commonly held view is that the response of the climate system to external forcing will tend to have a structure that is similar to the structure of one or a few of the system’s leading intrinsic modes of variability. Though intrinsic modes are often identified in terms of the large-scale upper-tropospheric circulation, typically they extend to the lower troposphere and influence the regional statistics of surface conditions and synoptic disturbances. Thus relating forced responses to intrinsic modes is a step toward understanding the regional character of that part of the forced response that directly affects society. The application of this idea that is perhaps most often cited is the reaction of the climate system to increasing anthropogenic greenhouse gas emissions, with the similarity in circulation trends during the end of the twentieth century to the North Atlantic Oscillation (Hurrell 1995) being taken as an example of a forced trend that resembles a prevalent natural mode. Interestingly, though there is broad agreement that one should expect a correspondence between forced and intrinsic structures, there is marked disagreement as to why one would expect the climate system to have this property. In this paper we report on an analysis that has been carried out of a climate change numerical experiment whose design facilitates an evaluation of the merits of competing explanations.

Though there may be many possible reasons that the forced response and leading intrinsic modes of the climate system should have similar structure, there are two explanations that are most often cited and we will restrict our analysis to them. The first is based on the observation that analysis of the linearized governing equations indicates, when considering system anomalies, there will be patterns whose natural frequencies are comparatively small and thus easier to force with quasi-steady sources than all other structures. Because of their low frequencies and ease of stimulation these patterns will be prevalent in analyses of both intrinsic and forced variability. Branstator (1985b, 1990) and Newman et al. (1997) are examples of investigations that have discussed this concept. Given the successful application of linear theory to many problems in interannual variability, including the formation of observed intrinsic patterns of atmospheric variability (Simmons et al. 1983) and the response to external forcing (Hoskins and Karoly 1981; Branstator 1985a), the assumption that linear reasoning is sufficient for understanding any similarity between forced and intrinsic anomalies is well grounded. The known importance of nonlinear flux terms in maintaining some prominent intrinsic patterns (Lau 1988; Branstator 1992) may seem to call into question the adequacy of this linear explanation, but it is known that these terms can be well approximated linearly (Branstator and Haupt 1998; Peng et al. 2003; Jin et al. 2006), so this basic explanation for similarity is viable even for such patterns.

The second prevalent explanation for the similarity in structure of forced and intrinsic variability relies on mechanisms that cannot be represented by linear equations. Rather this explanation asserts that nonlinearities in the system equations produce preferred states that show up as features of enhanced probability in probability density functions (PDFs) so that these PDFs are non-Gaussian. Variability in such a system largely consists of transitions between the preferred states, and the structure of the variability will be the difference in the two states involved. An especially pronounced example of this occurs when nonlinearities lead to multiple equilibria (Charney and DeVore 1979) and associated multiple local maxima in PDFs (Molteni 1996). According to this second explanation, when the dynamics of a system are under the influence of preferred states, its reaction to an external source can come about as a redistribution of the likelihood that the system will be in each of its prominent states while the structure of those states does not change (Ghil 1987; Palmer 1999). This means the structure of the response will be a weighted combination of pairwise differences between preferred states. So, again, there will be a match in the structure of the intrinsic and forced variability.

That such a paradigm needs to be taken into account to understand similarities between forced trends and modes of variability is supported by studies that have found evidence of non-Gaussianity (Kimoto and Ghil 1993a; Smyth et al. 1999; Hannachi 2007) or even multiple local maxima in PDFs of atmospheric states (Corti et al. 1999). There is some question about the statistical significance of these results (Hsu and Zwiers 2001; Stephenson et al. 2004), but Berner and Branstator (2007) found that in extremely long simulations with an atmospheric general circulation model (GCM) statistically significant non-Gaussianity in planetary wave states was present. Their distributions had only one local maximum but were consistent with the GCM having two overlapping clusters of states and so could not be used to rule out the possible relevance of the second explanation for the similarity between intrinsic modes of variability and forced responses.

In our investigation we have examined a climate change experiment in which the structure of circulation trends is clearly related to prominent structures that occur as a result of intrinsic variability. Our main aim has been to determine whether either of the two reasons outlined above effectively explains this correspondence. Note that, even though both of the explanations for the correspondence between forced and intrinsic structures use the term “modes of variability,” they use it in different ways. In the first explanation a mode of variability refers to a commonly occurring anomalous *pattern* with no intrinsic amplitude. We will refer to this as the *linear/pattern explanation*. In the second explanation the term refers to the transition between two preferred *states*. These states have a particular amplitude and the structure of the mode of variability is the structure of the difference between these states. We will refer to this as the *nonlinear/state explanation*. One can think of the main aim of our study as being to determine whether either usage of the term “modes of variability” is appropriate for describing the behavior of the climate change experiment we have examined.

The experiment that we have considered is the climate change experiment reported by Selten et al. (2004; hereafter SBDK). As described in section 2, this is an ensemble numerical experiment produced by a coupled ocean–atmosphere–land–sea ice model to simulate the reaction of the climate system to projected increases in greenhouse gases and aerosols. SBDK reported that the ensemble mean surface temperature trend during that part of the experiment that simulates the years 1940–2000 closely matches the observed temperature trend, and some ensemble members match the observed evolution of the North Atlantic Oscillation (NAO) during this period. Furthermore, SBDK found that the circulation trend in this experiment has a distinctive structure that is similar to the Circumglobal Waveguide Pattern (CWP). Branstator (2002) found the CWP was an intrinsic mode of nature and of some models, and Brandefelt and Koernich (2008) have reported that in climate change experiments the response of many models has this same structure. All of these factors indicate this experiment is a good candidate for studying the relationship between forced climate trends and intrinsic modes. What makes the SBDK experiment especially well suited for our study is that it is composed of an unusually large ensemble of realizations all reacting to the same external forcing. Because of this experimental design we have been able to take a systematic sequence of steps to evaluate the merits of the explanations outlined above for satisfactorily explaining why forced and intrinsic modes have similar structure, at least in this particular experiment. This has involved determining whether both the intrinsic variability and forced response in the experiment are consistent with each of the explanations.

The first three steps that we have taken address the question of whether trends in the SBDK experiment are consistent with the linear/pattern explanation. 1) As described in section 2, we have documented trends in the ensemble of model solutions for various fields of interest, including those that describe the atmospheric circulation, synoptic eddy activity, and surface weather conditions. 2) We have used empirical orthogonal function (EOF) analysis to find prominent patterns of intrinsic interannual variability, and then we have employed regression to find variations in circulation, synoptic eddies, and weather associated with these patterns. By comparing these to the forced trends we have quantified to what degree the climate change is shaped by the prominent intrinsic variability. These results are discussed in section 3. 3) We have investigated whether linear atmospheric dynamical mechanisms are sufficient to produce the dominant intrinsic patterns by finding the leading patterns generated by stochastic excitation of the linearized governing equations. These results, presented in section 4, together with additional experiments with the atmospheric component of the coupled model used by SBDK, indicate that the linear/pattern explanation matches many aspects of the simulated trends.

To complete our analysis we have taken additional steps to determine whether the nonlinear/state explanation is equally successful or perhaps explains features that are not consistent with linear behavior. 4) We have calculated and analyzed PDFs of projections onto the leading patterns of intrinsic variability. Then we have determined whether they have features that are consistent with the nonlinear/state explanation both in terms of their configuration before the effects of increased greenhouse gas forcing are felt (section 5) and in terms of how they change in reaction to enhanced greenhouse gases (section 6). We have been especially interested in whether there are indications of preferred *states*, whose probability, rather than circulation structure, changes as the climate changes. 5) Interestingly, though there is some suggestion of a redistribution between preferred states, the most important feature that is inconsistent with linearity is that the structure of the system’s preferred patterns and states changes as the climate changes. Therefore we have striven to understand what has produced this change in modal structure (section 7). In section 8 we summarize our results and conclude that the SBDK experiment is an example of a system in which both linear and nonlinearities influences can be seen in the interplay between modes of variability and forced trends, but the linear/pattern explanation appears to be adequate for understanding most prominent aspects of the similarity in intrinsic and forced structures.

## 2. Experimental climate change trends

The numerical experiment we have examined (Amman et al. 2003; SBDK) is a standard business-as-usual climate change experiment produced by the National Center for Atmospheric Research’s (NCAR's) Community Climate System Model (Boville et al. 2001), version 1.4 (CCSM1.4; Otto-Bliesner and Brady 2001). CCSM1.4 is a coupled climate model, including a T31 18 sigma level atmospheric component and a 3° oceanic component, that has been used in numerous studies. In this experiment it was used to simulate the years 1940 through 2080 by forcing the years 1940 through 1999 with time-dependent greenhouse gas concentrations, sulfate and volcanic aerosols distributions, and solar radiation as they were observed to occur during those years and by forcing year 2000 and beyond with sulfate aerosols and solar radiation fixed at their year 2000 values and with greenhouse gases continuing to increase at a rate similar to that in the Special Report on Emissions Scenarios (SRES) A1b scenario.

The primary reason for producing a simulation like the SBDK experiment is to learn how the earth’s climate might change in the future. But this results in a dilemma. For, by definition, climate properties involve averaging functions of state over a long enough time period to produce a statistically stable quantity. But when dealing with a changing climate, this means that many properties of interest are not well defined since they are changing during the averaging period. As Leith (1978) has explained, the solution to this dilemma is to consider an ensemble of systems, each experiencing the identical external conditions, and to average across the ensemble in addition to averaging in time. Thus shorter time intervals are needed to produce stable statistics, and if enough realizations are available the time interval can be short enough for the forcing (and thus the climate) to be considered to be stationary. Following this approach, SBDK repeated their experiment 62 times, with each realization being identical except for the initial condition. Because of this we have been able to determine many statistical aspects of the model climate that cannot be carefully considered in experiments that employ smaller ensembles. This is particularly pertinent for the evaluation of the statistical behavior of quantities whose intrinsic time scale is long.

As SBDK have reported, the secular changes in forcing in this experiment produce distinctive changes in the model climate. Given our goal of seeking a dynamical interpretation of the regional aspects of these changes, we have examined the manifestation of these changes in various fields of dynamical and meteorological interest. In each case we have used the ensemble- and time-average field during years 2051–80 of the experiment minus the ensemble- and time-average field during years 1941–70 to represent the climate change. We refer to these time intervals as the *final* and *initial periods* of the experiment, respectively. Anticipating that the structure and dynamics of climate change are likely to be a function of season, we have calculated these differences for a single season, namely December–January–February (DJF). We will refer to such differences as the model *trend*.

The trends calculated in this way and discussed below have all been based on averages of all 62 realizations, but we have also calculated the trend in each field of interest two additional times. These additional estimates were arrived at by dividing the 62 realizations into two halves. Except where we state otherwise, every prominent feature in trends and other statistical properties calculated from one-half were also present in trends and statistics derived from the other half, thus establishing that climate change properties presented in the following paragraphs result from forcing variations and not from internal processes.

We have found that from a dynamical perspective one of the most revealing trend fields is the seasonal mean 300-hPa streamfunction. As shown in Fig. 1a, the anomalies induced in this field have a prominent zonal wavenumber 5 component at the latitude of the subtropical jet. A field that emphasizes these jet-trapped anomalies even more is the north–south wind component, whose trend at 300 hPa is displayed in Fig. 1b. Though the zonal wavenumber 5 scale of trends is most notable in the upper troposphere, this structure is also detectable nearer the surface, as can be seen in Fig. 1c’s plot of trends in 850-hPa streamfunction. SBDK showed that there are also distinctive trends in this experiment in various fields in the Southern Hemisphere, but with our focus on relating trends to modes of variability, we have chosen to only show the Northern Hemisphere portion of the trends. The leading modes we have found and that are discussed later are confined to the Northern Hemisphere.

Using the same approach we have also examined trends in quantities that have a more direct bearing on the statistics of surface weather. One of these is near-surface temperature, which we have represented by temperature at 850 hPa. The trend in this field is displayed in Fig. 1d, which shows the strongest trends are concentrated in northern polar latitudes and over the eastern portions of the Northern Hemisphere continents. This distribution is similar to the cold ocean–warm land (COWL) pattern that Wallace et al. (1996) found was associated with recent long-term variations in temperature, though the land maxima in their pattern are shifted toward the centers of continents. Another weather indicator is synoptic eddy variability, which we have found by calculating the variance of geopotential height at *σ* = .324, Fourier filtered to retain periods between 1 and 7 days. As seen in Fig. 1e, Northern Hemisphere anomalies in this quantity are largest over the ocean basins, where the storm tracks are strongest. The final weather-related field we have examined is precipitation, shown in Fig. 1f. It reflects shifts in precipitation associated with both midlatitude changes in storm-track position and tropical changes in convection.

## 3. Intrinsic patterns of variability

As mentioned in the introduction, the experimental trend that SBDK found in the 300-hPa streamfunction and meridional wind fields (our Figs. 1a,b) is reminiscent of an intrinsic pattern of variability that Branstator (2002) discovered in atmospheric general circulation model integrations and in nature. He called it the Circumglobal Waveguide Pattern because it appears to result from the waveguiding effect of the mean subtropical jet on planetary waves with zonal wavenumbers greater than or equal to 5, an effect that has been discussed by Branstator (1983) and Hoskins and Ambrizzi (1993). Branstator (2002) determined that the meridional wind in the upper troposphere was an especially good field to use when isolating this intrinsic phenomenon.

Given this resemblance we have checked whether the CWP is an intrinsic pattern of DJF averages in the SBDK experiment. So that the statistics of the model climate are essentially stationary in this part of our analysis, we have restricted our attention to the initial period of the experiment. Furthermore, so that only internal variations are considered, we have subtracted the ensemble average of each winter from each of the 62 samples for that winter. Then, as a means of isolating the leading *patterns* of intrinsic variability, we have found the leading EOFs of the resulting 62 × 30 charts. When we have done this for 300-hPa meridional wind, the two leading eigenvectors that resulted are the patterns displayed in Fig. 2. These patterns make it clear that zonal wavenumber 5 disturbances confined to the subtropical jet are prominent components of internal variability. The two vEOFs are in spatial quadrature, so together they span a family of wave 5 patterns, but as Branstator (2002) found with other AGCMs, the longitudinal phasing of vEOF1 (which explains almost 3 times as much variance as vEOF2) is by far the most commonly occurring orientation of members of this family. Since vEOF1 is very similar in terms of structure and variance dominance to the pattern that Branstator identified, we consider it to be the model’s CWP.

Comparing the structure of the model CWP with the climate change trend in Fig. 1b, we see that not only do the zonal scales and latitudinal placement of the centers of action match but the longitudinal phasing also agrees. Hence the regional character of the trend, at least in this field, appears to be strongly influenced by the presence of this intrinsic pattern. To be more quantitative about how much of the trend in the meridional wind can be explained by the intrinsic CWP, we have found the degree to which it by itself can approximate the Fig. 1b trend. To do this we first used regression to find the distribution of winds that are naturally associated with vEOF1. This entailed projecting our collection of 62 × 30 charts of internal DJF 300-hPa meridional wind anomalies onto vEOF1 and then regressing this time series against the gridpoint values of these same wind fields to form a map of dimensional meridional winds. We did this in such a way as to produce the meridional winds associated with a one standard deviation vEOF1 event. Next, we scaled the resulting map by *c*, which is the projection of the Fig. 1b trend onto vEOF1 divided by the standard deviation of vEOF1 projections.^{1} The result of this calculation is shown in Fig. 3b. Clearly this single pattern is a very good approximation to the trend. This conclusion is confirmed by the fact that the field in Fig. 3b explains 80% of the spatial variance north of 15°N in Fig. 1b.

Branstator (2002) found that temporal variations in the polarity and amplitude of the CWP are associated with variations in other quantities. For example, though the global character of the CWP is most pronounced in the upper troposphere, its amplitude covaries with winds throughout the depth of the troposphere in some regions. Furthermore one might expect it to be associated with statistical anomalies in synoptic activity given the close linkage between large-scale circulation anomalies and bandpass eddy statistics that has been described in many papers (Lau 1988; Branstator 1995; Whitaker and Sardeshmukh 1998). Considering the trends in lower-tropospheric and weather statistics noted in the previous section, it is natural to wonder whether the regional structure of those trends corresponds to the structure in these quantities expected from the influence of the CWP. To investigate this possibility we have calculated anomalies associated with natural variations in the CWP in the same fields whose trends are shown in Fig. 1.

We have found anomalous conditions associated with internal variations in the CWP during the initial period of the experiment by using the same linear regression approach used to construct Fig. 3b. For each field of interest *f* we have isolated the internal variations for DJF averages during each simulated year by subtracting the ensemble mean for that year from each realization. Then for each gridpoint we have calculated the regression coefficient that relates internal variations of *f* at that location to a one standard deviation intrinsic variation in the model CWP. Finally we have scaled these regression maps by *c* so that they correspond to the distribution of *f* associated with the trend in vEOF1.

Figure 3a shows the 300-hPa streamfunction produced by this sequence of operations. Considering the close physical relationship between streamfunction and meridional wind, and considering the similarity between the CWP and the trend in meridional wind, it is not surprising that the 300-hPa streamfunction anomalies associated with the CWP are very similar to the trends in this field (Fig. 1a). Indeed the CWP explains 71% of variance in Fig. 1a north of 15°N. A similar correspondence between CWP-related variability and trends follows when we repeat this analysis but for the circulation in the lower troposphere. As seen in Fig. 3c, the 850-hPa streamfunction anomalies associated with intrinsic variations in the CWP are primarily composed of two features in the northern ocean basins that represent equivalent barotropic extensions of lows in the upper troposphere. These same vertical extensions are seen in the trends (Fig. 1c). By contrast, the regional structure of near-surface temperature trends are not nearly as strongly influenced by the CWP, as can be seen in the intrinsic anomalies in 850-hPa temperature associated with the CWP (Fig. 3d). They have much less similarity to the trends (Fig. 1d) than do the other fields we have examined, though some trend features, including positive anomalies near Hudson Bay and northwest Africa, may be a result of the CWP. It is also possible that the cooling caused by the CWP over northwestern Eurasia is the reason that the model warming trend over Asia is farther east than in the observed COWL pattern.

Just as impressive as the above influence of the CWP on trends in mean fields is its influence on the statistics of synoptic-scale variability. A comparison of Fig. 3e, which displays the intrinsic anomalies in upper-tropospheric bandpass geopotential variance associated with the CWP, and model trends in this quantity (Fig. 1e) indicates that most of the prominent trends in the storm tracks are controlled by the CWP. Both the generally southward shift in eddy variance and regional features like the dipoles just west of the Greenwich meridian and just east of the date line appear to owe their existence to the influence of the CWP. North of 15°N, the CWP-associated anomalies explain almost 40% of the structure of synoptic eddy variance trends. A second indicator of synoptic activity is midlatitude rainfall. Its trends are also influenced by the CWP. This is seen in Fig. 3f, which shows anomalous precipitation during intrinsic CWP events. The main features (on the eastern sides of the ocean basins and southwest of Mexico) also appear in the trends (Fig. 1f). Interestingly, there are also some strong precipitation anomalies in the tropics associated with intrinsic CWP variability, but without further analysis it is not apparent whether these should be thought of as causing or being caused by the CWP.

## 4. Modes from stochastically driven linear dynamics

The results of the previous section strongly suggest that regional features in trends in many fields in the SBDK climate change experiment are related to similar features in a prominent intrinsic pattern of the model climate. But these results to do not help distinguish the dynamical processes that produce that intrinsic pattern nor do they rule out the possibility that that *pattern* is actually a reflection of preferred *states*. To address these issues we have tested whether simple linear dynamics by themselves are sufficient to understand why the CWP is the most prominent pattern of intrinsic variability in the model used by SBDK.

Using the linearized barotropic vorticity equation, Branstator (2002) showed that the waveguiding effect of the climatological subtropical jet could lead to structures like vEOF1 being trapped in the jet. But his experiments did not establish that these structures should necessarily be more prominent than other structures, and it is not known how sensitive his results might be to the specific basic state used in the computation. To find out if these linear mechanisms do select wavenumber 5 for prominence in the wintertime subtropical jet of the system used by SBDK, we have taken the adiabatic equations from the atmospheric component of CCSM1.4 and linearized them about the time-averaged, ensemble-averaged DJF state from years 1941–70 of the SBDK experiment. To solve these equations, the same discretization methods used by Branstator (1990) were employed with a horizontal truncation of R15 and 9 vertical levels. To make up for the essential stabilizing role of the nonlinearities that are missing in the linear system, we have added Rayleigh and Newtonian damping terms with characteristic time scales of 1 day in the bottom layer, 2 days in the next layer, and 7 days above that. Then, as a means of finding the prominent intrinsic features of this linear system, we have used a method similar to that used by Branstator (1990). This involved finding many steady solutions of the model, with each forced by a different steady source of vorticity whose three-dimensional structure was determined by assigning gridpoint values based on draws from a Gaussian distribution of random numbers. Before being employed as forcing, the zonal average was removed from each forcing.^{2} Since the forcing has no organization in this experiment, any organization found in the solutions must be from the linear dynamics and the effect of the basic state.

The stationary wavenumber arguments of Hoskins and Karoly (1981) and the result of Held (1983) suggest that the zonal mean structure of the wintertime climate may be sufficient to produce the jet trapping and scale selection that characterize the CWP, so we first considered solutions linearized about the zonal mean component of the SBDK climate. Figure 4 shows the first two EOFs of the meridional wind fields at *σ* = 0.336 that come out of 1000 randomly forced solutions. Clearly linear dynamics do lead to zonal wavenumber 5 disturbances in the vicinity of the jet being important intrinsic quasi-stationary perturbations.

Though the previous calculation produced patterns with many key attributes of the leading patterns of the SBDK experiment, the properties of the linear model in that calculation mean that the leading two patterns it produces cannot have two characteristics that are seen in the SBDK experiment. One of the two patterns cannot have more variance than the other and their meridional structures cannot be a function of longitude. To see whether these characteristics may result from the action of zonal asymmetries in the SBDK climate, we have repeated the stochastic linear calculation but with the full three-dimensional climate mean for the background field.

Figure 5 shows vEOF1 and vEOF2 for this calculation. The basic zonal wavenumber 5 structure remains but now the two additional properties just noted are produced by the linear mechanisms. Both patterns meander farther north over the North Pacific and over western Europe than at other longitudes just as happens in CCSM1.4. And now one pattern explains more variance than the other, with the longitudinal phasing of the dominant pattern being the same in the linear and SBDK experiments. Note, for example, that for both CCSM1.4 and the linear model, vEOF1 has a lobe over the northern Red Sea while that lobe is north of the Arabian Sea in both versions of vEOF2. To be sure, there are some features of the linear vEOFs that do not agree with the CCSM1.4 behavior. For example vEOF1 is not as dominant and its longitudinal dependence is not as prominent as it should be in the linear calculation. One possible reason for these discrepancies is the lack of feedbacks from transient eddy fluxes in this calculation. We have not tried to incorporate these feedbacks into our linear experiments, because we believe the above results offer compelling evidence that the linear/pattern explanation can account for most important features in the intrinsic CWP of the SBDK experiment.

The stochastic linear calculations also appear to be consistent with the regional structure of the trend in the SBDK experiment, but they offer no insight into the sign of trend anomalies. To complete our investigation of whether the trend is consistent with the linear/pattern interpretation of the similarity between forced response and modes of variability, we have considered whether there is any means by which one particular sign of the leading intrinsic mode may be selected. The mechanism we have tested concerns whether the external forcing trends in the experiment may have induced a heating trend that stimulates the leading linear pattern. As explained by Branstator (1985b) it is not necessary for such a source to project more strongly onto the leading mode than onto other modes for the leading mode to dominate the response. But there must be a source that has some projection onto that mode.

Examination of Fig. 1 suggests that the source that may be forcing the CWP trend is heating associated with the prominent trend in precipitation in the far western tropical Pacific (Fig. 1f). This is the strongest anomaly in the tropics, and the lower-tropospheric streamfunction trends (Fig. 1c) have the signature of a wave train emanating from near this location, entering the waveguide and exciting the CWP. Figure 6 displays a test of this hypothesis in which we have placed a steady heat source in the linear model configured as it was for the stochastic tests with a three-dimensional basic state. This heat source approximates the western Pacific precipitation trend by being concentrated in the midtroposphere and decreasing linearly from its maximum value of 5 K day^{−1} at 10°N, 150°E to a value of zero on an ellipse with a zonally oriented semimajor axis of 2000 km and a semiminor axis of 1000 km. The meridional wind response at *σ* = .336 displayed in Fig. 6 shows that such a source does produce a response that projects strongly on vEOF1 from the stochastic experiment (Fig. 5a). Though structurally correct throughout the Northern Hemisphere this response is markedly weak eastward of 90°W. Again, this discrepancy may result from the lack of transient eddy feedbacks. Other calculations we have done but are not shown here indicate they help to pump energy into the CWP over the North Atlantic.

The linear response of Fig. 6 does indicate that heating near the western Pacific precipitation trend would tend to force the CWP and give it the sign of the trend in the SBDK experiment. Further support for this possibility is given by another experiment we have performed in which we have used the atmospheric component of CCSM1.4 and forced it with a sea surface temperature (SST) anomaly near this same location. This experiment consisted of a 9-yr control integration in which this model was forced by climatological values of SST and sea ice taken from years 1961–90 of the SBDK experiment and a 9-yr anomalous experiment in which the model forcing was altered by the addition of a Gaussian shaped SST anomaly centered at 2°N, 148°E. The central value of this anomaly was 3 K, and it had a standard width of 15°. It induced a precipitation anomaly of placement and structure much like the SST anomaly. As Fig. 7 indicates, the anomalous 300-hPa meridional wind anomalies induced by this heating have many of the important features of the trend and intrinsic CWP we have been discussing.

## 5. Preferred states

Though the results of the previous two sections generally support the linear/pattern explanation for the similarity of intrinsic and forced anomalies in the SBDK experiment, we have also considered whether preferred states and associated transitions, as hypothesized by the nonlinear/state explanation, may be responsible for, or at least influence some aspects of, this similarity. To determine whether there are such states we have examined phase space distributions of intrinsic variability in the SBDK experiment to see whether there are multiple regions of high density and whether there are indications of the non-Gaussianity that should be associated with such regions.

Even though we have more data than are typically available for studying distributions of low-frequency states, we have been restricted to considering distributions in two dimensions so that our results are robust. Even then, we have needed to employ analysis methods that use either substantial smoothing or parametric fits so that we could meet a certain standard of robustness. That standard is the same one we used for investigating trends; namely, that only features found in analyses of both of two halves of our dataset are considered robust.

The basic analysis method we have used to investigate whether there are signatures of non-Gaussianity in the distribution of states in the SBDK experiment is to estimate and examine PDFs of its states. We have done this using a kernel method (Silverman 1986; Kimoto and Ghil 1993a) that employs a bivariate Gaussian kernel with a window width *h* = 0.2. (With *h* = 0.1 our standard of robustness is not satisfied.)

We have focused our attention on the plane defined by projections onto vEOF1 and vEOF2 (Fig. 2). This plane has the advantage of including the CWP, which we found to be largely responsible for the experiment trend we have been investigating. Furthermore, by construction, intrinsic anomalies should be stronger in this plane than any other, a property that should enhance the possibility that nonlinearities will influence their behavior.

Projecting intrinsic variations in DJF averages of 300-hPa meridional winds for the initial period of the experiment onto vEOF1 and vEOF2, normalizing by the standard deviation of the vEOF1 projections, and then calculating the kernel-based PDF for these data, we have produced the distribution displayed in Fig. 8a. The marked departures of contours in this plot from ellipses and the displacement of the maximum from zero are indications that this distribution is not Gaussian. A measure of this non-Gaussianity is

which is the relative entropy of the kernel-based PDF *ρ* relative to the reference distribution *ρ _{G}*, which is the bivariate Gaussian that best fits the data. Kleeman (2002) has argued relative entropy is a useful measure of how two distributions differ from each other. We find

*R*= 0.037.

If the PDF had multiple modes, that is, multiple local maxima, this would be a strong indication of the appropriateness of the nonlinear/state usage of “modes of variability.” The PDF does not have multiple modes, but this does not exclude the possibility of multiple modes in some higher dimension that strongly overlap when projected onto our plane. Or, as Majda et al. (2006) have pointed out, a non-Gaussian PDF without multiple peaks can be consistent with a system that has multiple metastable states that lead to overlapping clusters of high probability states. As shown in Itoh and Kimoto (1996), these states could be related to remnants of true multiple attractors of the system when its parameters are modified. To check whether possibilities like these have any merit we have seen how well the PDF can be approximated as a mixture of a small number of bivariate Gaussians. This was accomplished using the expectation-maximization (EM) algorithm described by Berner and Branstator (2007). EM algorithms (Titterington et al. 1985; Haines and Hannachi 1995; Smyth et al. 1999) are designed to maximize the log likelihood of a mixture of Gaussians for a particular dataset.

When we fit data in the vEOF1–vEOF2 plane with a mixture of just two Gaussians, we found that it produced a distribution that is a very good approximation to the kernel-based PDF. This can be seen in Fig. 8b, which shows the individual components as well as the distribution that results when they are combined. The position of the distribution maximum, the pronounced skewness along a line through this maximum and the origin, and a secondary ridge of enhanced probability in the first quadrant are all reproduced by the mixture. An objective measure of the accuracy of this approximation is that the log likelihood of the two-Gaussian mixture is −13 428 while the log likelihood of the best fit single Gaussian is −13 563.

When we tried approximating the data in the vEOF1–vEOF2 plane with more than three Gaussians, our robustness standard was not met; component parameters varied markedly for different subsamples. On the other hand, a three-component mixture could not be completely ruled out as a possible approximation to the PDF. Positions and shapes of the individual Gaussian components bore clear similarities in the two subsamples that we considered but varied too much to make the tests of sensitivity of component parameters to climate change reported in the next section viable. For this reason, as well as the fact that the log likelihood of the three-component mixture only improved to −13 390, we have decided to restrict our attention to the two-Gaussian mixture as a model of the distribution of states in the vEOF1–vEOF2 plane that captures its most dominant non-Gaussian features. Note that Berner and Branstator (2007) reached a similar decision based on a more comprehensive analysis of internal variability in a different AGCM. And Hannachi (2007) has come to the conclusion that a two-Gaussian mixture is the appropriate approximation to daily reanalysis charts of Northern Hemisphere winter 500-hPa geopotential heights.

It is possible that even more distinctive indications of non-Gaussianity are present in other directions. To check this possibility, we considered other planes. So that we only considered planes that include the trend from the SBDK experiment, each of the planes we considered included a normalized version of Fig. 1b as one direction. The other direction was formed from linear combinations of the five leading EOFs of intrinsic 300-hPa DJF meridional wind variability. We calculated 1000 such patterns by assigning random coefficients drawn from a centered uniform distribution to each of the five vEOFs. These patterns were then modified in such a way as to be orthogonal to the trend pattern. We then calculated *R* for kernel-based PDFs for each of the corresponding 1000 planes. The pattern that, together with the trend pattern, defines the plane with the largest relative entropy is shown in Fig. 9a. Except for an east–west-oriented dipole over western North American, the structure of this pattern is very similar to vEOF2. (In our robustness test this dipole was present in only one subsample, while the wave 5 anomalies in the subtropics were present in both subsamples.) Considering that the trend pattern is very similar to vEOF1, it is not surprising then that the PDF for this plane (Fig. 9b) is quite similar to Fig. 8a. This similarity was further verified when we fit the data in this plane with a two-component Gaussian (not shown). The position and orientation of these Gaussians is very similar to those in Fig. 8b. The relative entropy for the kernel-based PDF in this plane is 0.041. Thus we have found that there are planes with distributions that are more non-Gaussian than the vEOF1–vEOF2 plane, but even the most non-Gaussian of these is not very different, nor much more non-Gaussian, than this plane. Hence we think it is reasonable to view the two Gaussian components of the vEOF1–vEOF2 plane as being representative of two (mildly) preferred states of the system.

## 6. Changes in state distribution during climate change

The mode of variability implied by the two preferred states derived from the mixture analysis is the difference of the component centroids and thus has a very prominent vEOF1 component, as can be inferred from Fig. 8b. This finding is consistent with the nonlinear/state explanation for the similarity between modes of variability and forced response. To further consider the applicability of that explanation, we have examined how the PDF of the previous section changes as the climate changes in that experiment. Figure 10a shows the kernel-based PDF for the final period of the SBDK experiment, where the axes pertain to the same vEOFs used in analysis of the initial period. Interestingly, although the mean of the distribution shifts in the positive vEOF1 direction, consistent with the trend projecting onto positive vEOF1, the new distribution does not appear to be a simple shift in the distribution. For example, the region of frequently occurring high-amplitude anomalies in the second quadrant of the initial PDF (Fig. 8a) does not shift with the mean but rather straddles the −vEOF1 axis in the final PDF. Changes like these are summarized in Fig. 10b, which is a plot of the final PDF minus the initial PDF.

To see whether the mixture model from the previous section has any bearing on the changes portrayed in Fig. 10, we have fit a two-Gaussian mixture to data from the final period. As shown in Fig. 11a, such a mixture is a good approximation to the final PDF. And if we subtract the initial two-Gaussian mixture (Fig. 8b) from this distribution, the result (Fig. 11b) captures the essence of the change in the PDF (Fig. 10b). Hence if we can understand the change in the two-Gaussian mixtures we can understand the change in the complete PDF.

Comparing the pairs of Gaussian components in Figs. 11a, 8b, we see that during the final period there has been an enhancement in the population of the right component. As we have noted, the nonlinear/state paradigm assumes a system with multiple modes in its PDF will react to an external forcing through a redistribution of probability among the modes without changing (the phase space position of) the modes. Taking into account that the two-component mixture we use to approximate the initial period of the SBDK experiment has the form

where the *G*s are bivariate Gaussians and the *α*s are weights, if such a description is valid for the SBDK experiment, then the final PDF should be well approximated by (1) with different *α*s but the same *G*s. We have checked this possibility by finding those *α*s that best fit the two-Gaussian mixture from the final period of the experiment. Here best is measured by relative entropy. This optimal mixture is shown in Fig. 11c. Note that by construction the position and shape of the two components are identical to those in Fig. 8b. The error in this distribution (Fig. 11d), as given by subtracting the actual final two-Gaussian mixture from this optimized mixture, is somewhat smaller than the error one would make by assuming there had been no change in the distribution (Fig. 11b), but the improvement is modest. By contrast, when we have simply shifted the initial mixture by the actual change in the mean values of the distributions, the resulting distribution (Fig. 11e) has much smaller errors (Fig. 11f). This last approximation amounts to assuming that the effect of the trend on variability is purely additive.

Though the linear shift appears to more accurately approximate the change in distribution of states than the redistribution implied by the nonlinear/state paradigm does, its errors are substantial. This failing of linearity should not be surprising given some of the attributes we noted in comparing the initial PDF to the final PDF. Further examination of these PDFs and their mixture components suggests an alternative to the two mechanisms tested in Fig. 11. We have already pointed out that the second quadrant ridge of enhanced probability in the initial PDF rotates counterclockwise to a new position in the final PDF. Similarly, if one compares the centroids of the Gaussian components for the initial and final distributions (given by closed and open circles, respectively, in Fig. 11a), they also experience a rotation from the climate change. To quantify this rotation we have taken the initial PDF, shifted it by the change in the mean that occurs between the beginning and final state, and then rotated it about this new mean through that angle that produces a distribution that most closely matches the final PDF. Again we used a relative entropy measure of closeness. This optimally shifted and rotated PDF is displayed in Fig. 12 together with its error, which is now substantially smaller than it was for either of the standard mechanisms.

The position and orientation of the original coordinate axes when they are shifted and rotated in the same fashion that the PDF has been transformed in Fig. 12a are plotted in that figure as dashed lines. Comparing the transformed abscissa to the dashed line in Fig. 10a, which gives the direction in the final PDF with the most variance, we see that to a large extent the optimal rotation has simply aligned the direction of maximum variance in the initial PDF with the direction of maximum variance in the final distribution. Said another way, the rotation of the initial PDF that transforms it into a distribution that nearly matches the final PDF represents the fact that the pattern with maximum variance has changed as the climate changed. Consistent with this observation, when we have plotted vEOF1 for the final 30 yr of the SBDK experiment (not shown) we have found that its structure is much like its structure during the initial 30 yr but shifted slightly eastward. This shift is easy to see when we plot in Fig. 13e vEOF1 for the final period minus vEOF1 from the initial period. This eastward shift means that the leading pattern has taken on a contribution from the original vEOF2, as is consistent with the sense of the rotation of axes in Fig. 12a.

## 7. Mechanism producing change in patterns of variability

To understand why the patterns of variability have changed, we have hypothesized that the same linear mechanism we found to be largely responsible for the initial leading patterns could also be responsible for the final leading patterns. That is we have hypothesized that the trend in the coupled model’s mean state changes the background state in such a way that the leading patterns produced by linear dynamics are modified and those modifications match the vEOF changes seen in the experiment.

Tests with the same linear model used in section 4 appear to confirm this hypothesis. When we have linearized it about the three-dimensional average state from the final period of the SBDK experiment and forced it with the same random functions employed in section 4, the solutions have a leading vEOF for meridional winds at *σ* = .336 that is much like linear stochastic vEOF1 for the initial basic state (Fig. 5a) but shifted eastward (not shown). As with vEOF1 from the final stage of the SBDK experiment, this shift can be highlighted by plotting in Fig. 13f the new stochastic vEOF1 minus stochastic vEOF1 for the initial basic state. Comparing Figs. 13e,f we see that the mechanism represented in the linear stochastic experiment is sufficient to reproduce most features of the change in structure of the leading pattern.

As further confirmation of our hypothesis, we have also considered two intermediate stages of the SBDK experiment. Figures 13a,c display how the structure of vEOF1 for two 30-yr periods, roughly one-third and two-thirds of the way through the experiment, differs from its initial structure. Taken together with Fig. 13e, they show that vEOF1 gradually shifts more and more eastward as the experiment progresses. When we have repeated the linear stochastic excitation experiments but for basic states based on these intermediate stages, we have found that the leading linear pattern also gradually shifts eastward (Figs. 13b,d,f).

## 8. Conclusions and discussion

We have found that for the climate change experiment reported by SBDK, the frequently expressed idea that the structure of forced regional climate trends are likely to be influenced by the structure of intrinsic modes of variability is valid. Not only are mean circulation trends in the upper troposphere, where the modes are often defined, affected by these modes, but so are trends in lower-tropospheric mean circulation and weather statistics. As to whether this is because of either of two explanations that are prominent in the literature, or perhaps because of some third explanation, our results give mixed conclusions.

Our results clearly show that the linear/pattern explanation, that involves prominent *patterns* produced by mechanisms that can be represented linearly, effectively explains much of the behavior seen in the climate change experiment. The climate trend has the same structure as the Circumglobal Waveguide Pattern, which is a pattern that appears as a leading structure in a linear-based analysis of intrinsic variability. That same pattern turns out to be a prominent pattern of steady-forced variability in the linearized equations of motion, presumably because its natural frequency is near zero. There is a trend in the tropical heating in the climate change experiment that is capable of forcing this linear mode and giving it the sign it has in the trend. PDFs of state vectors that include the trend pattern as a component are quasi-Gaussian and have a single maximum. To a first approximation the change of the PDF as the climate changes is a simple additive shift. All of these attributes are consistent with a linear/pattern framework for understanding the correspondence between modes of variability and regional trend structure.

As successful as the linear/pattern interpretation is, there are aspects of the climate change experiment that are inconsistent with this framework. PDFs for some phase space planes that include the trend pattern are not purely Gaussian. Rather they are well approximated by a mixture of two overlapping Gaussians, a configuration that could be indicative of there being two preferred system *states*. Moreover, the changing structure of the PDF as the climate changes is not completely describable as a simple shift.

Interestingly, these inconsistencies with the linear/pattern explanation are not entirely consistent with the most commonly stated hypothesis for the similarity of intrinsic and forced structures that involves preferred states resulting from nonlinearities. This hypothesis, which we called the nonlinear/state explanation, involves transitions between nonchanging preferred states and assumes that the response to external forcing consists of a redistribution of probability between these states. Though our results did show a suggestion of such a redistribution as the climate reacts to the greenhouse gases, on closer examination the departure from the linear/pattern framework appeared to be better described as a gradual change in the structure of the leading intrinsic patterns of the system. A defining characteristic of this changing structure is that the direction of maximum variance rotates in phase space planes. Hannachi (2007) has noted a rotation in PDFs during the reanalysis era, but he analyzed different patterns and a different field, namely 500-hPa geopotential heights. Brandefelt (2006) has also reported evidence of changing prominent patterns in a climate change experiment with a coupled model, but she was not able to successfully find a dynamical interpretation for this change.

Our results suggest that the change in pattern structure in the SBDK experiment is a reaction of the linear modes to the modifications to the climatological mean state that are forced by the ramp up in greenhouse gases. The linear stochastic results, together with the similarity between intrinsic linear modes and circulation trends reported in section 3, indicate that the experiment provides an example of a two-way interaction between intrinsic patterns and trends. The intrinsic (linear) patterns of the undisturbed climate shape the structure of the mean reaction to the change in climate forcing and that change in the mean produces modifications in the intrinsic patterns. Note that even though we used linear analysis and a linear model to decipher these interactions, they do constitute nonlinearities. Further analysis suggests that one more nonlinear interaction between patterns and trends may be at work. From Fig. 10a we see that the mean of the distribution of states at the end of the climate change experiment falls closer to the direction of maximum variance for this part of the experiment than to the abscissa, which is the direction of maximum variance for the beginning of the experiment. Thus it appears that the trend structure may be influenced by the modified leading pattern, a pattern it helped to produce.

Finally, it is worth putting the nonlinear influences, no matter how they are described, in perspective by noting that their impact is comparatively weak; it is only because of the large dataset that we have studied that we were able to distinguish the influence of the nonlinearities from the dominant linear characteristics of intrinsic variability and climate response. Using the terminology of the introduction, this means that for our study the use of “modes of variability” that pertains to *patterns* is more appropriate than the usage that refers to transitions between preferred *states*. On the other hand, though it did not appear to be particularly important for understanding climate response in this experiment, the possible importance of non-Gaussianity and preferred system states cannot be ruled out for other situations. For example, our results do not exclude the possibility that, in cases where the behavior of regional rather than global structures is more important than in the SBDK experiment (e.g., Kimoto and Ghil 1993b), or where system behavior on shorter time scales is of interest (e.g., Dole and Gordon 1983), or where a system with more prominent preferred states is considered (e.g., Palmer 1999), the nonlinear/state mechanism may be more relevant.

## Acknowledgments

This investigation began while the first author was a visiting scientist at KNMI and continued while the second author was a visiting scientist at NCAR. The study was further supported by Contract NA04OAR4310061 from the NOAA Climate Program Office program in Climate Variability and Predictability, and by the Office of Science (BER), U.S. Department of Energy, Cooperative Agreement DE-FC02-97ER62402. The figures were generated by A. Mai, and the text benefited from comments offered by three anonymous reviewers.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

## Footnotes

*Corresponding author address:* Grant Branstator, National Center for Atmospheric Research, P.O. Box 3000, Boulder, CO 80307-3000. Email: branst@ucar.edu

^{1}

Alternatively, we could have simply scaled vEOF1 by the projection of the trend in 300-hPa meridional wind onto vEOF1. We used the method described in the text because it is easily modified to find fields other than 300-hPa meridional wind associated with vEOF1.

^{2}

We have also used a time-marching method and forced the system with time-dependent random Gaussian forcing that was white in time. This gives similar results if one analyzes 30-day or longer averaged perturbations.