## Abstract

A nonlinear principal component analysis (NLPCA) is applied to monthly mean zonal wind observations from January 1956 through December 2007 taken at seven pressure levels between 10 and 70 hPa in the stratosphere near the equator to represent the well-known quasi-biennial oscillation (QBO) and investigate its variability and structure. The NLPCA is conducted using a simplified two–hidden layer feed-forward neural network that alleviates the problems of nonuniqueness of solutions and data overfitting that plague nonlinear techniques of principal component analysis. The QBO is used as a test bed for the new compact model of NLPCA.

The two nonlinear principal components of the dataset of the equatorial stratospheric zonal winds, determined by the compact NLPCA, offer a clear picture of the QBO. In particular, their structure shows that the QBO phase consists of a predominant 28.3-month cycle that is modulated by an 11-yr cycle as well as by longer cycles. The differences in wind variability between westerly and easterly regimes and between Northern Hemisphere winter and summer seasons and the tendency for a seasonal synchronization of the QBO phases are well captured.

## 1. Introduction

The familiar quasi-biennial oscillation (QBO) dominates variations of the zonal winds in the stratosphere (from about 16- to 50-km altitudes) around the equator. In general, the zonal wind at each pressure level alternates between easterly and westerly with a period varying between about 20 and 36 months with a relatively rapid transition over about 2 to 4 months. The wind regimes consistently propagate downward through the equatorial stratosphere. The wind amplitude, associated with the QBO, of about 20 m s^{−1} is nearly constant between 10 and 40 hPa and then reduces to about 5 m s^{−1} at 70 hPa. The overall wind oscillations show noticeable asymmetries: (i) the transition of easterly to westerly is more rapid than that of westerly to easterly, and (ii) the associated westerly shear zone, where westerly winds increase with height, descends more regularly and rapidly than the easterly shear zone, which sometimes stalls for several months between about 30 and 50 hPa (Naujokat 1986; Marquardt and Naujokat 1997).

The equatorial stratospheric zonal winds exhibit a rather complicated height–time structure. There has been considerable work to construct low-dimensional statistical fits that describe the key features of the overall oscillations of the zonal winds covering a depth of the equatorial stratosphere. For instance, the principal component analysis, also called empirical orthogonal function analysis, has been applied to the monthly mean zonal winds in the 10–100-hPa range of the equatorial stratosphere. It was found that the leading two principal components varied out of phase through the QBO period, and the leading two modes could describe the basic height–time structure of the zonal winds reasonably well (Wallace et al. 1993). A singular spectrum analysis, also called extended empirical orthogonal function analysis, was applied to those zonal winds as well (Fraedrich et al. 1993; Wang et al. 1995). Although the linear approaches provided a reasonably good representation of the QBO structure, some of the more interesting characteristic features of the wind oscillation are lost. For instance, the reconstructed wind series using the leading two modes from the principal component analysis tend to look much more sinusoidal in time than the actual wind oscillations, the shear zones are considerably smoothed out, and the asymmetries between easterly and westerly regimes are almost completely missed.

Recently, a three–hidden layer feed-forward neural network with a circular bottleneck node (Kirby and Miranda 1996) for nonlinear principal component analysis (NLPCA), called the NLPCA.cir model, was used to represent the monthly mean zonal winds observed between 10 and 70 hPa in the equatorial stratosphere (Hamilton and Hsieh 2002). Taking advantage of the periodic nature of the QBO, in this approach the observations are mapped to a single phase variable that advances between −*π* and *π* over each QBO cycle. This one-dimensional NLPCA.cir model solution fits the data more accurately than even the two-mode reconstruction from the principal component analysis. More importantly, this nonlinear, one-dimensional representation captures the qualitative features of the typical QBO cycle much better than the linear, two-mode reconstruction of the zonal winds.

One important benefit of a low-dimensional statistical fit to the equatorial stratospheric zonal winds is that a QBO phase can be determined that accounts for the overall variations of the zonal winds through the depth of the equatorial stratosphere. In the linear fit by the leading two modes from the principal component analysis, a time series of QBO phase is defined as the arctangent of the leading two principal components (Wallace et al. 1993). In contrast, the single variable in the NLPCA.cir model solution can be directly interpreted as the phase of the QBO (Hamilton and Hsieh 2002).

An estimate of QBO phase as a function of time is needed for at least two types of investigations: (i) characterizing the connections of the tropical stratospheric QBO with other aspects of the atmospheric circulations and (ii) studying the relation between the QBO and the annual cycle. The first of these research areas has a long history extending back at least to the work of Holton and Tan (1980), who found that in the northern extratropics the stratospheric circulation differed on average between easterly and westerly phases of the equatorial stratospheric QBO. The phase of the QBO has also been correlated with ozone variations and many aspects of meteorological variability in the stratosphere and troposphere (Baldwin et al. 2001). The possible relation between the QBO and the annual cycle has been discussed by various investigators (Dunkerton and Delisi 1985; Dunkerton 1990; Wallace et al. 1993). Interestingly, there have been a number of recent papers reexamining both of these aspects of the QBO. The connection between the phase of the tropical stratospheric QBO and the weather of the northern troposphere, an issue of possible importance for seasonal forecasting, has been examined based on observations (Thompson et al. 2002; Cai 2003). The problem of the annual synchronization of the QBO, in observations and simple mechanistic model experiments, has been revisited (Hampson and Haynes 2004). Almost all of these studies have used a definition of the phase of the QBO based on the time series of the zonal wind at a single level, which is often chosen to be 50 or 40 hPa. An exception is the examination of the annual synchronization of the QBO phase, which is determined from the leading two principal components of the linear fit (Wallace et al. 1993). By contrast, the NLPCA.cir model incorporates the full vertical structure of the observed zonal winds into the time series of a single variable interpretable as the phase of the QBO cycle (Hamilton and Hsieh 2002). As noted above, this nonlinear approach works quite well in many respects and better than the linear principal component analysis, but the restriction to essentially a constant amplitude in representing the wind oscillations is a significant limitation.

The constant amplitude produced by the NLPCA.cir model is imposed by the architecture of the model. This neural network model consists of five nonlinear and linear feed-forward layers of neurons. In sequence, they are input, encoding, bottleneck, decoding, and output layers, with three hidden layers sandwiched between the input and output layers. The NLPCA.cir model is given only two bottleneck neurons and constrained to project the two bottleneck series onto a unit circle. For this reason, the pair of bottleneck neurons is called a circular bottleneck node. Consequently, the projection allows only the variations in phase progression to be carried into the decoding neurons (at the expense of amplitude variations, which are lost). Only angular phase progression is extracted from the pair of transformed bottleneck series; their amplitude variations are ignored.

A straightforward approach to produce a nonlinear representation with variable amplitudes is simply to relax the constraint of a circular bottleneck node. This would effectively replace the NLPCA.cir model with the standard three–hidden layer feed-forward neural network for NLPCA originally proposed by Kramer (1991). However, it is now known that applications of the standard three–hidden layer model are subject to the technical problems of generating nonunique solutions and overfitting the data (Hsieh 2004). In this study, we adopt a newly developed, more robust approach of NLPCA using a simplified two–hidden layer feed-forward neural network, called the compact NLPCA model (Lu and Pandolfo 2008, manuscript submitted to *Neural Networks*, hereafter LP) (see Fig. 1) to represent and characterize the QBO. This model alleviates the problems of nonuniqueness of solutions and data overfitting that plague nonlinear techniques of principal component analysis.

In this paper, the compact NLPCA model will be applied to the monthly mean zonal winds observed at seven pressure levels in the equatorial stratosphere. The paper is organized as follows: Section 2 briefly describes the data under analysis. Section 3 introduces the compact NLPCA model. Section 4 selects an optimal two-component representation of the data. Section 5 analyzes the variations of the QBO phase defined by the two nonlinear principal components. Long period variations in the phase progression of the QBO are also investigated. Section 6 discusses the variations of the wind variability with phase and season. Section 7 shows the seasonal synchronization of the QBO phases and the variations of the phase speeds with phase and time, and section 8 offers conclusions.

## 2. Data

The data to be analyzed are the well-known monthly mean series of the stratospheric zonal winds (Naujokat 1986; Marquardt and Naujokat 1997), calculated from the twice-per-day measurements by balloons above Canton Island (2.8°N, 171.7°W) from January 1956 to August 1967, Gan (Island) in the Maldives (0.7°S, 73.1°E) from September 1967 to December 1975, and Singapore (1.4°N, 103.9°E) from January 1976 to December 2007. The operational balloon soundings are usually capped at 10 hPa, although the QBO appears at least up to 40-km altitude (Hamilton 1981). The monthly mean series at 10, 15, 20, 30, 40, 50, and 70 hPa (from about 30- to 20-km altitude) provided by Free University of Berlin (available online at http://www.geo.fu-berlin.de/en/met/ag/strat/products/qbo/index.html) are used in this study and are hereafter called the QBO winds for convenience. As described above, the QBO winds consist of seven time series, each having 52 yr of monthly data.

## 3. The compact NLPCA model

Our NLPCA is built on a neural network model. Hence, the terminology used to describe the technique follows the language of neural networks. A neural network (Fig. 1) is composed of a number of layers, with each layer consisting of an ensemble of neurons. The neuron is the basic unit of a neural network. The type of data (e.g., input data) assigned to the neurons of a layer identify that layer (e.g., the input layer). In general, a neural network at least has an input layer of neurons for data to enter the neural network for analysis and an output layer of neurons for the model to produce its final analysis output. This output is the approximation of the original data built by the neural network model during its evaluation of the nonlinear principal components of the data. The layers sandwiched between the input and output layers are called hidden layers, and their neurons are correspondingly called hidden neurons. It is within these hidden layers that the evaluation of the nonlinear principal components takes place. LP describe how to construct a compact NLPCA model with a reduced number of hidden layers that transforms classical NLPCA into a quasi-objective multivariate statistical analysis technique.

The compact NLPCA model (Fig. 1) has four layers with two hidden layers of neurons. In sequence, they are input, bottleneck, decoding, and output layers. Each variable of the dataset under study is represented by a specific input neuron. The data are called input signals {*p _{in}*}, where

*p*denotes the signal from input neuron

_{in}*i*at sampling point

*n*,

*i*= 1, … ,

*I*, with

*I*being the number of input neurons, and

*n*= 1, … ,

*N*, with

*N*being the number of samples. The number of input neurons

*I*equals the number of variables.

The compact NLPCA model built to analyze the QBO winds will have seven input neurons (*I* = 7), with each input neuron corresponding to the zonal wind time series at a specific pressure level. The parameter *n* refers to a sample in time and *N* (= 52 yr × 12 months yr^{−1} = 624 months) is the number of data in each of the seven monthly mean series.

To produce the NLPCA of the QBO winds, the input data must flow through the neural network in the following manner. First, all input signals are combined through a nonlinear transfer function to generate the time series *u _{kn}* associated with each bottleneck neuron of the first hidden layer:

where *u _{kn}* denotes the signal from bottleneck neuron

*k*at sampling point

*n*;

*A*is the weight of bottleneck neuron

_{ki}*k*to signals from input neuron

*i*; and

*k*= 1, … ,

*K*, with

*K*being the number of bottleneck neurons. The hyperbolic tangent function has become the usual choice of nonlinear transfer function for neural networks and it is also used in this study. Dimension reduction is achieved when signals pass through the bottleneck neurons. Therefore,

*K*, the number of bottleneck neurons, is always less than that of input neurons (

*K*<

*I*). The time series from each bottleneck neuron is taken as a nonlinear principal component. For example, given two bottleneck neurons (

*K*= 2), the model will extract two nonlinear principal components from the data under analysis. (In the following sections,

*u*= {

*u*|

_{n}*n*= 1, … ,

*N*} denotes one and

*υ*= {

*υ*|

_{n}*n*= 1, … ,

*N*} denotes the other of the two nonlinear principal components of the QBO winds.) The bottleneck weight

*A*could be considered as a scaling factor in projecting signals from input neuron

_{ki}*i*onto axis

*k*of a

*K*-dimensional plane. In general,

*A*varies with both

_{ki}*k*and

*i*.

Second, all bottleneck signals are combined again through a nonlinear transfer function to generate the time series *y _{mn}* associated with each decoding neuron of the second hidden layer:

where *y _{mn}* denotes the signal from decoding neuron

*m*at sampling point

*n*;

*C*is the weight of decoding neuron

_{mk}*m*to signals from bottleneck neuron

*k*;

*c*is the bias of decoding neuron

_{m}*m*, and

*m*= 1, … ,

*M*, with

*M*being the number of decoding neurons, which is adjustable for an optimal fit of output to input signals. Usually

*C*varies with both

_{mk}*m*and

*k*. The signals {

*y*} are called decoding signals because the equations for

_{mn}*y*are inverse nonlinear transforms of the bottleneck signals; that is, they effectively decode the nonlinear principal components. The purpose of this decoding becomes clear after the last transfer to the output layer.

_{mn}Finally, all decoding signals are forwarded to each output neuron through a linear transfer:

where *q _{in}* denotes the signal from output neuron

*i*at sampling point

*n*and

*D*is the weight of output neuron

_{im}*i*to signals from decoding neuron

*m*. The number of output neurons is equal to the number of input time series. The output signals {

*q*} are the reconstructions of the input time series based on the representation of the input by the nonlinear principal components estimated in the bottleneck layer of the model. This representation is made optimal through an iterative procedure that minimizes the difference between the output signals and the input time series. In this iterative procedure, the proper number of decoding neurons is seldom known beforehand; in practice, a series of values of

_{in}*M*is tried out for NLPCA to find the optimal fit between output and input signals.

In our case, the compact NLPCA model should have seven output neurons when applied to the QBO winds. The signal series from each of the output neurons is a low-dimensional approximation of that from the corresponding input neuron. Effectively, the input–bottleneck part of the neural network projects the input time series nonlinearly to a low-dimensional space (the space of the nonlinear principal components) and the bottleneck–decoding–output part performs an inverse transform projection back to the original space represented by the model output. Because signals are always fed forward from one layer of neurons to the next, the neural network is said to be feed forward. It is also described as auto-associative because output signals are fitting to input signals and thereby output neurons are in pairs with input neurons. The neural network will hereafter be referred to as an *I*-*K*-*M*-*I* model for the specific architecture.

If the neural network is given only one decoding neuron, all the output series will be linearly related to each other because of the linear transfer from decoding to output signals. If all transfer functions are set to be linear, a feed-forward neural network will degenerate to perform a linear principal component analysis (Bourlard and Kamp 1988; Baldi and Hornik 1989). To focus on nonlinear analysis, only the neural networks with certain nonlinear transfer functions and two or more decoding neurons (*M* ≥ 2) are used in this study.

In NLPCA, appropriate values of all the weight and bias parameters discussed above are determined by a minimization of the mean square error of output to input signals. The minimization is accomplished by means of a hybrid procedure described in the appendix.

## 4. Nonlinear representation of the QBO winds

The standard deviations of the seven time series of the QBO winds vary from 6.4 to 20.0 m s^{−1}. Proper scaling of the QBO winds helps to attain an optimum performance of the NLPCA. Accordingly, the QBO winds are normalized by removing their 52-yr means and dividing them by their overall standard deviation (i.e., the standard deviation of the total data) and then a 7-2-*M*-7 model is applied, where *M* is varied from 4 to 20 with 40 runs for each *M*. Two bottleneck neurons are assigned to the neural network to extract two nonlinear principal components from the QBO winds and allow the variability in both phase progress and amplitude value of the QBO winds from cycle to cycle to be effectively captured.

The model output, called simulation, which also consists of seven time series, is multiplied by the overall standard deviation. The 52-yr means are then added to the corresponding times series and the root-mean-square error is calculated for comparison with the original data and shown in Fig. 2a. The lowest root-mean-square error of the 40 runs for each *M* is 3.51, 3.29, 3.16, 3.08, 3.01, 2.95, 2.91, 2.89, 2.87, 2.84, 2.84, 2.81, 2.81, 2.79, 2.78, 2.78, and 2.78 m s^{−1} for *M* = 4, … , 20, respectively. It reduces relatively quickly from 3.51 to 2.81 m s^{−1} as *M* is varied from 4 to 15 and then remains between 2.81 and 2.78 m s^{−1} when 15 ≤ *M* ≤ 20 (Fig. 2b). The 7-2-15-7 model, containing the fewest parameters among the models offering the lowest root-mean-square errors, is optimal in terms of least squares error and fewest model parameters.

The time series of the 7-2-15-7 model simulation (Fig. 3, solid line) closely approximate the corresponding time series of the QBO winds (Fig. 3, dots), with correlations 0.974 ≤ *ρ* ≤ 0.994 for pressure levels between 10 and 50 hPa (and *ρ* = 0.866 for 70 hPa, where the original data are rather noisy).

The standard deviation *σ* of the errors of the 7-2-15-7 model simulation to the QBO winds varies between 1.99 and 3.60 m s^{−1} among the seven pressure levels (Fig. 4). Compared to the normal distribution with expectation zero and variance equal to that of the errors of the 7-2-15-7 model simulation (Fig. 4, curve), the frequency distribution of the errors of the 7-2-15-7 model simulation has a higher kurtosis and little skewness (Fig. 4, vertical lines). Because the 7-2-15-7 model simulation is a very satisfactory fit to the QBO winds, the two 7-2-15-7 nonlinear principal components are an effective two-component representation of those winds.

For comparison, a standard linear principal component analysis was also performed. The wind series reconstructed from the leading two linear modes (not shown) was characterized by correlations 0.929 ≤ *ρ* ≤ 0.983 for pressure levels between 10 and 50 hPa and *ρ* = 0.806 for 70 hPa. The standard deviation of the errors of the linear two-mode reconstruction to the QBO winds varies between 3.65 and 5.93 m s^{−1} among the seven pressure levels (Fig. 5). Compared to a normal distribution with expectation zero and variance equal to that of the errors of the two-mode reconstruction (Fig. 5, curve), the frequency distribution of the errors has a lower kurtosis and/or a noticeable skewness (Fig. 5, vertical lines). The nonlinear two-component representation of the QBO winds is significantly better than the linear two-mode reconstruction.

## 5. Variations of the QBO phase

Figures 6a and 6b show the time series of the two 7-2-15-7 nonlinear principal components, *u* and *υ*, for the full 52 yr. The time series of each component undergoes a QBO oscillation, of course. The amplitude of the oscillation varies only modestly from period to period, but the period shows more obvious cycle to cycle changes. This can be seen also in Fig. 6c, which shows a measure of the phase defined by Ψ = arctan(*υ*/*u*) (−*π* ≤ Ψ ≤ *π*). The variability of QBO period from cycle to cycle has been reported in earlier studies, but these have generally based their measure of the QBO phase on the zonal wind at a single level (Quiroz 1981; Dunkerton and Delisi 1985; Naujokat 1986; Maruyama and Tsuneoka 1988).

### a. Mean frequency of the QBO

The long-term mean period of the QBO has been estimated in earlier studies. The fairly recent paper of Baldwin et al. (2001) uses 4 decades of data and finds the mean period is slightly longer than 28 months. Here we address this issue using the results of our nonlinear principal component analysis.

The periodic phase Ψ (Fig. 6c) shows a strong tendency for the phase to increase linearly with time through each cycle. So an accumulating phase Φ, defined as the progression of the phase Ψ but with each new cycle starting at (2*n*_{cycle} − 3)*π*, with *n*_{cycle} varying from 1 to 24, is calculated. The first cycle of Φ starts from −*π* because the phase Ψ advances from −*π* to *π* in each cycle. The accumulating phase Φ progresses predominantly in a linear fashion (Fig. 7a), characterized by the linear regression Φ* _{L}* = 2

*πn*

_{month}/28.3 + 0.2

*π*, where

*n*

_{month}is the number of months from the initial sampling point and starts from 1. It indicates that the phase of the QBO winds advances predominantly at a speed of 2

*π*/28.3 months. We performed a simple Fourier analysis of the time series of the winds at each pressure level and of the winds simulated with the NLPCA. The magnitudes of the Fourier components are displayed in Fig. 8, where the strong power at 22 × 1/52 = 0.423 cpy (equivalent to a 28.4-month period) is very evident.

It is noticeable that the phase analysis can produce an accurate phase speed that matches the periodicity of the QBO winds obtained by their Fourier analysis. Furthermore, one can get this phase speed through the phase analysis even based on datasets that have different lengths (e.g., *N* = 480 or *N* = 510). The coincidence between the values of the QBO period determined through our phase analysis and the Fourier transformation of the wind data could be fortuitous. This is because the Fourier transformation of the QBO winds with 52 years of data has a frequency resolution of 1/52 cpy. The frequency of 22 × 1/52 = 0.423 cpy, which is equivalent to 28.4-month period, is followed by a frequency of 23 × 1/52 = 0.442 cpy, equivalent to 27.1-month period. Hence, contrary to our phase analysis, the Fourier transformation analysis of the zonal wind time series cannot differentiate between 28-month and 28.3-month periods.

### b. Frequency modulation of the QBO phase progression

Recent studies have raised the possibility of systematic quasi-decadal modulation of the QBO phase progression. Salby and Callaghan (2000) showed a tendency for the westerly wind regimes at 45 hPa to be shorter near solar maxima over the 1956–96 era and suggested a connection of the QBO phase with cycles of sunspot number of about 11 yr. Hamilton (2002) noted that the apparent connection was rather weak for the solar maximum near 1981 and was less apparent when data from the early 1950s and late 1990s were included. Hamilton suggested that some mechanism unrelated to solar activity might be responsible for generating a quasi-decadal variation in the QBO phase progression. Recently, Fischer and Tung (2008) analyzed the equatorial zonal wind observed from 1953 to 2007 at each of the seven pressure levels between 10 and 70 hPa to reexamine the modulation of the 11-yr sunspot cycle on the QBO period. They found that the QBO period changed little with height and that although the QBO period seemed generally lengthened (shortened) in solar minima (maxima) from the 1950s to the 1980s, an opposite or in-phase relationship started in the 1990s.

Here we examine these issues using the phase defined on the basis of the NLPCA. Figure 7b shows the time series of a phase residue {*ϕ _{n}*}, defined as the actual phase in any month minus that obtained by assuming the mean phase progression for the full period. To emphasize the 11-yr cycle and its persistent appearance in the QBO phase, 97 subsets are each consecutively sampled from the phase residue for 44 years (i.e.,

*ϕ*, … ,

_{t}*ϕ*

_{t+527}for each subset, where

*t*= 1, … , 97) and Fourier analyzed in a straightforward way. Figure 7c shows the magnitude of the Fourier components (densely overlapping dots at each frequency; only the amplitude coefficients at frequencies between 0 and 1 cpy are shown because those at higher frequencies are much smaller). The spectra of Fourier amplitude coefficients of these subsets closely resemble each other, indicating that the spectral pattern is a robust feature of the phase residue. In the spectra of Fourier amplitude coefficients, the phase residue has two adjoining, prominent peaks at 0.023 cpy (= 44-yr period) and 0.045 cpy (= 22-yr period), with another prominent peak at 0.091 cpy (= 11-yr period) (Fig. 7c). Nonetheless, because we have 52 years of data, the spectral peaks at 0.023 and 0.045 cpy represent modulations with periods close to the length of the data. The strength of these oscillations might not have been sampled properly. Unfortunately, the signal of the multidecadal period is not well documented. In an ordinary meteorological record spanning a few decades, this signal can hardly be picked out from strong short-term signals, and it is easily mixed into a nontrivial long-term trend. Therefore, a discussion of these possible modulations cannot be done in a meaningful way. The following discussion will focus on the dominant peak at 0.091 cpy that corresponds to a significant modulation of the QBO phase at a period of 11 yr.

Similar to the result of Hamilton (2002), the time series of the harmonic component of 0.091 cpy (Fig. 7c; asterisk) follows the ups and downs of the phase residue during the 1960s and 1970s but completely misses the variations in phase residue from the mid-1980s to the present (Fig. 7b; line). The 11-yr cycle by itself is not enough to describe the overall modulation of the QBO phase.

An interaction between two oscillations at frequencies *f* and *g* can produce another two oscillations, one at frequency *f* − *g* and the other at frequency *f* + *g*. The interaction can be expressed mathematically as

where *a* and *b* are constant offset parameters, *f* and *g* denote frequencies, and *t* is time. The interaction between a biennial oscillation at *f* = 0.50 cpy and a QBO at approximately *g* = 0.41 cpy could produce an oscillation at *f* − *g* = 0.09 cpy or an 11-yr period. The signal of the biennial oscillation would disappear when *b* is negligible. The phase residue obtained from the 7-2-15-7 model solution might be such a case.

However, in the tropical regions, radiative ozone heating can introduce an 11-yr sunspot cycle to stratospheric dynamics. Nonetheless, the two possible origins of the 11-yr cycle, one from the variation of solar irradiance and the other through the interaction of relevant oscillations in the atmosphere, do not exclude each other. The possibility exists that both mechanisms contribute to the modulation of the QBO.

### c. Statistical significance of the QBO phase features

The essential linear increase with time of the accumulating phase Φ and the dominant decadal fluctuations of the phase residue *ϕ* are robust features present in the two nonlinear principal components. These features all appear and are almost the same among the solutions of the 7-2-*M*-7 model runs with close results of the root-mean-square error. For instance, 40 runs were conducted for each of the 6 models 7-2-{15, 16, 17, 18, 19, 20}-7, that is, 240 runs in total. Their resultant root-mean-square errors all fall between 2.78 and 2.88 m s^{−1} (Fig. 2a). The accumulating phases defined by the two nonlinear principal components of each of the 240 individual solutions closely resemble each other (not shown), having correlations larger than 0.99. These 240 accumulating phases have a mean of 2*π*/28.3 months and a standard deviation of 2*π*/0.0018 months for the slope coefficients of their linear regressions. Compared with the value of the mean, the standard deviation is negligible. The corresponding phase residues (not shown) are correlated at correlations between each other greater than 0.95 among 97% of the 240 series. From the point of view of reproducibility, the predominant linear increase at a speed of 2*π*/28.3 months and the dominant decadal modulations at 11-, 22-, and 44-yr cycles of the phase residue are robust and key features of the QBO phase.

The Fourier amplitude coefficients of the phase residues at frequencies higher than 6 × 1/44 ≈ 0.14 cpy are significantly smaller than those at 1/44 ≈ 0.023, 2 × 1/44 ≈ 0.045, and 4 × 1/44 ≈ 0.091 cpy, as seen in the spectra of Fourier amplitude coefficients of the phase residue from the 7-2-15-7 model solution (Fig. 7c).

## 6. Variance in the variability of the QBO winds

To describe variations of the QBO winds, a reference state needs to be constructed. This reference state will be a general structure of the QBO without perturbation, in which the amplitude of oscillations of the zonal wind at each pressure level is constant at each phase. This general structure should be extracted starting from the two 7-2-15-7 nonlinear principal components because they represent the QBO winds. This is realized by using a 2-2-*M ^{ψ}*-2 model to generalize the two 7-2-15-7 nonlinear principal components; that is, the two 7-2-15-7 nonlinear principal components are the input and target series of the 2-2-

*M*-2 model, where

^{ψ}*M*is varied from 2 to 9 with 10 runs for each

*M*. The superscript

*ψ*denotes that the signal series entering the decoding neurons are those of the sine and cosine of the angle series formed by the arc tangent of the two bottleneck series of the 2-2-

*M*-2 model, rather than the bottleneck series themselves. Hence, in this model, the two signal series entering the decoding neurons have a constant amplitude. Subsequently, the output signals of the 2-2-

^{ψ}*M*-2 model are constant at each phase of the angles, as shown in the next paragraph.

^{ψ}The lowest root-mean-square error of the 10 runs for each 2-2-*M ^{ψ}*-2 model is 0.074, 0.066, 0.064, 0.063, 0.061, 0.061, 0.060, and 0.060 for

*M*= 2, 3, … , 9, respectively. Except for the relatively quick reduction of the root-mean-square error as

*M*is varied from 2 to 3, reduction of the error with further increase of

*M*is rather slow. The significant change at

*M*= 3 in the rate of decrease with increasing

*M*of the lowest root-mean-square error among the 10 model runs and the relatively small number of parameters contained in the model with

*M*= 3 are logical conditions for the choice of

*M*. So the 2-2-3

*-2 model solution is optimal based on the combined effect of least squares error and fewest model parameters. Presented in a scatterplot, the generalized nonlinear principal components [i.e., the two series of the 2-2-3*

^{ψ}*-2 model simulation (Fig. 9a; densely overlapping crosses)] pass through the 7-2-15-7 nonlinear principal components (Fig. 9a, dots) as a smooth and closed curve. Entering the generalized nonlinear principal components into the decoding neurons of the optimal 7-2-15-7 model, the corresponding seven output series form a general height–phase structure of the QBO winds (Fig. 9b).*

^{ψ}Relative to the generalized nonlinear principal components, the two 7-2-15-7 nonlinear principal components appear more variable in one half of their phase cycle, from −0.75*π* counterclockwise to 0.25*π*, than in the other half from 0.25*π* to −0.75*π* (or equivalently from 0.25*π* to 1.25*π*). Referred to the general height–phase structure, the phase range from −0.75*π* to 0.25*π* corresponds to a duration when westerlies prevail through the pressure levels and easterlies are absent or weak except for their buildup at upper levels around −0.25*π*. In contrast, the phase range from 0.25*π* to 1.25*π* corresponds to a duration when easterlies dominate through the pressure levels and westerlies are absent or weak except for their buildup at the upper levels around 0.75*π*. The variability in the amplitudes of the nonlinear principal components shows marked differences between the two phase regimes of the QBO winds. The variability is greater during the westerly phase than the easterly phase.

The variance also changes with season. After being stratified according to seasons, the data points in the scatterplot of the two 7-2-15-7 nonlinear principal components diverge frequently from the closed curve of the generalized nonlinear principal components in the half phase cycle from −0.75*π* to 0.25*π* in February–April (Fig. 10a) and November–January (Fig. 10d). This indicates that the nonlinear principal components and hence the QBO winds are rather variable in this phase range during the boreal winter. In contrast, the data points in May–October (Figs. 10b,c) tend to be close to the closed curve of the generalized nonlinear principal components through the full phase cycle. The nonlinear principal components and hence the QBO winds are quite stable during the austral winter whether the winds are blowing west or east.

These changes in wind variability with phases and seasons are consistent with the annual cycle of Rossby waves in the stratosphere and the asymmetry between northern and southern atmosphere in the interaction of the equatorial stratospheric winds and the linear planetary Rossby waves. Consider the general height–phase structure of the QBO winds (Fig. 9b). Whereas an easterly or westerly starts to reach its full strength at 10 hPa, a wind of opposite direction is strong at 40 hPa. Hence, composites of a field variable based on the equatorial zonal wind at 10 hPa will be nearly opposite to those based on the equatorial zonal wind at 40 hPa. The westerly regime from 70 to 40 hPa coincides with the QBO phases between −0.75*π* and 0.25*π* when westerlies are prevailing in the equatorial stratosphere and easterlies are absent or weak most of the time. From linear wave theory, the presence of equatorial stratospheric westerlies allows a considerable amount of planetary Rossby wave energy to penetrate into the tropics, leading to mitigated wave activities in the extratropics and intensified perturbations of the equatorial stratospheric winds, and hence greater variability there. In contrast, the easterly regime from 70 to 40 hPa corresponds to the QBO phases between 0.25*π* and 1.25*π*, when easterlies are dominating in the equatorial stratosphere and westerlies are absent or weak most of the time. This blocks a considerable amount of planetary Rossby wave energy from entering the tropics, leading to enhanced wave drag on the polar westerly vortex and leaving the equatorial stratospheric winds relatively undisturbed.

## 7. Seasonal synchronization of the QBO phases

### a. Early summer phase transitions

A recent analysis based on a record of the monthly mean zonal wind near the equator at 50 hPa alone showed that the easterly onset occurs more frequently in May–July and the westerly onset in May–June than the other months of the year (Baldwin et al. 2001). This annual synchronization of the phases of the QBO has been linked to seasonally varying tropical upwelling (Hampson and Haynes 2004). It will be shown that this phase synchronization is seen in the two 7-2-15-7 nonlinear principal components; we will also see how it varies throughout the year.

The data points of the two 7-2-15-7 nonlinear principal components appear, in the scatterplot, to cluster around some angular phases, and these angular positions change with the seasons of the year (Fig. 10). In more detail, the distribution of data points according to the periodic QBO phase and calendar month (Fig. 11) shows a tendency for a seasonal synchronization of the phases. If the data points were distributed evenly, each phase interval would contain six to seven events in a month (52 samples per month over eight phase intervals of 0.25*π*). In fact, each QBO phase has some particular, consecutive months in which it is more recurrent than in the other months of the year.

For instance, a QBO phase of −*π* occurred 8 or 9 times a month during January–April. The number reduced to 4 or less in the other months of the year. A QBO phase of −0.25*π* appeared 7 to 13 times a month from October to May of the next year, but only 3 or 4 times per month during the other months of the year. Moreover, the QBO phase tends to advance with the particular months. The QBO phases of −*π*, −0.75*π*, and −0.25*π* are most recurrent during January–April, May–October, and October–May, respectively. Then the QBO phases of 0, 0.25*π*, 0.5*π*, and 0.75*π* are most recurrent during February–June, June–August, July–October, and October–December, respectively. It is followed by the QBO phase of *π*, equivalent to −*π*, being most recurrent during January–April. The QBO phase spans 21 (April to December of the next year) to 29 (January to May of the third year) such months in a cycle.

In reference to the observation noted at the beginning of this section—that at 50 hPa, the onset of either the easterly or the westerly phase of the QBO is more likely to happen during the period of late spring to early summer—confirmation comes directly from the results of our NLPCA. Figure 9 indicates that transitions at 50 hPa occur around phases of −0.75*π* and 0.25*π*. Figure 11 shows that these phases recur more often in late spring and early summer. A similar result was obtained by Huesmann and Hitchman (2001), using composite maps of the two QBO phases built using their new QBO index. An advantage provided by NLPCA is, for instance, that the probability of a given transition occurring in a particular month can be calculated using Figs. 9 and 11 because NLPCA supplies a time series of the QBO phase. Analysis based on composite maps is limited in this respect because a composite map can only represent a temporal average.

### b. QBO phase speeds

If the QBO winds had the same phase speed all the time, the distribution of data points should have been no different between months and phases. Actually, the phase speed varies because the phase residue fluctuates with time (Fig. 7b). For comparison between phases as well as months, the individual phase speeds are calculated from the accumulating phase as ΔΦ* _{n}*/Δ

*n*

_{month}= (Φ

_{n+1}− Φ

_{n−1})/2 and stratified according to the periodic QBO phase and calendar month (Fig. 12).

All the phase speeds are between −0.03*π* and 0.27*π* per month. The phase speeds around the QBO phase of 0, when easterlies blow at the upper levels and westerlies exist at lower levels, and around the QBO phase of *π*, when westerlies blow at the upper levels and easterlies exist at the lower levels (Fig. 9b), are mostly below 0.1*π* per month through the year. In contrast, the phase speeds around the QBO phase of −0.5*π*, when westerlies dominate through the pressure levels, fall below 0.1*π* per month during July–August and then rise to 0.1*π* to 0.3*π* per month during September–June. The phase speeds around the QBO phase of 0.5*π*, when easterlies prevail through the pressure levels, fall below 0.1*π* per month during December–March, but with relatively few events in these months, and then rise to 0.1*π* to 0.3*π* per month during April–November.

Naturally, a QBO phase with a nearly zero speed in a month would persist in this month or a following month. A QBO phase with a high speed would have low recurrence because the phase will quickly advance to the next phase. These associations between the phase recurrence and phase speed are seen in Figs. 11 and 12. The phase speeds at the QBO phases of ±0.5*π* are usually higher than 0.1*π* per month (Fig. 12). Accordingly, these QBO phases rarely recurred in a month, except for those months near which the previous QBO phase had a relatively high recurrence (Fig. 11). In contrast, the QBO phases 0 and *π* were often nearly stationary (Fig. 12). They recurred more than other QBO phases and occurred preferentially during February–June and January–April, respectively (Fig. 11). In the general height–phase structure (Fig. 9b), a QBO phase of 0 corresponds to the westerly-to-easterly transition near 40 hPa and a QBO phase of *π* corresponds to the easterly-to-westerly transition near 40 hPa.

Hampson and Haynes (2004) showed that the seasonally varying tropical upwelling could tune the QBO phases toward an annual synchronization. However, the present analysis indicates that each of the QBO phases had some specific consecutive months during which it occurred preferentially compared to the other months of the year. The preferred months change with the QBO phase and all seasons have such months, regardless of whether the tropical upwelling is strong or weak. In addition, according to the phase speeds (Fig. 12) and the general height–phase structure (Fig. 9b), the processes that synchronize the QBO exert different effects on different QBO phases. For instance, in July–August, they significantly slow down the phase speeds if westerlies are prevailing through the pressure levels, but not if easterlies are dominating there. This suggests that waves, propagating upward from the troposphere, could be involved in regulating the phase speeds of the QBO winds. Some particular waves in July–August provide eastward momentum flux to maintain the westerlies in the equatorial stratosphere, but they thereby decelerate the easterlies.

### c. Uneven descent rate of QBO phases

Kinnersley and Pawson (1996) considered that the upwelling branch at the equator of the meridional Brewer–Dobson circulation could slow down the descent of the QBO easterly regime more than the QBO westerly regime. The Brewer–Dobson circulation is believed to be forced by breaking and dissipation of the planetary Rossby waves in the polar stratosphere (Holton et al. 1995; Randel et al. 2002; Hood and Soukharev 2003). The slow descent of easterly winds around 30 to 50 hPa was found to take place mostly during Northern Hemisphere winter when the dissipation is strongest (Holton et al. 1995). Similarly, the 7-2-15-7 model solution shows that the QBO phases of 0.5*π*, 0.75*π*, and *π*, when easterly winds dominate between 70 and 30 hPa (Fig. 9b), are more recurrent during October–December (e.g., event number = 18 in December) than other times of the year (e.g., event number = 8 in June; see Fig. 11). The QBO phases of −0.75*π*, −0.5*π*, and −0.25*π*, when westerly winds prevail between 70 and 30 hPa, are, as a group, distributed uniformly throughout the year (e.g., event number = 24 in December compared to 20 in June; see Fig. 11). Hence, the time series of the QBO phase, calculated from the nonlinear principal components, suggests that the upwelling branch of the Brewer–Dobson circulation slows down the descent rate of the easterlies in Northern Hemisphere winter, as suggested by Kinnersley and Pawson (1996).

## 8. Conclusions

The new compact NLPCA model (LP) is applied to the monthly mean zonal winds observed from January 1956 to December 2007 at seven pressure levels between 10 and 70 hPa in the stratosphere near the equator. A nonlinear, two-dimensional statistical fit of the dataset of the zonal winds is produced. The two nonlinear principal components very satisfactorily represent the seven time series and characterize the key features of the QBO more precisely than the first two linear principal components of the same dataset. A summary of the key features investigated in this study is presented below.

First, the two nonlinear principal components offer a phase series with a quasi-periodic cycle. The calculated accumulating phase increases with time essentially in a linear fashion at a rate of 2*π*/28.3 months, which is also the predominant oscillation frequency of the observed and simulated zonal winds. This rate is the essential phase speed of the zonal wind oscillation. In other words, 28.3 months is a fundamental period of the QBO. The phase residue, after the subtraction of the linear component from the accumulating phase, fluctuates with a predominant 11-yr period as well as longer periods. The existence of modulations by cycles of longer periods could be the reason why the modulation of the QBO phase cannot be completely explained by using the 11-yr sunspot cycle alone. Also, whether the 11-yr cycle mirrors sunspot activity is arguable; the cycle could originate from the interaction between the QBO and a biennial signal.

Second, the two nonlinear principal components show clearly that the zonal winds in the equatorial stratosphere are subject to perturbations, and therefore are more variable, when mean westerlies dominate through middle or low levels of the equatorial stratosphere in the boreal winter. The zonal winds are relatively undisturbed whenever mean easterlies dominate through the middle or low levels of the equatorial stratosphere as well as in the austral winter, no matter what direction the mean zonal winds are blowing. These relationships between the variability and direction of the zonal winds in the equatorial stratosphere are consistent with the annual appearance and disappearance of the planetary Rossby waves in the northern extratropical stratosphere, the modulation of the waveguide by the equatorial zonal winds of the QBO, and the asymmetry in the strength of the planetary Rossby waves between the Northern and Southern Hemispheres.

Finally, the QBO phase defined by the two nonlinear principal components exhibits a tendency for seasonal synchronization. In general, a QBO phase recurs more frequently during some consecutive months than the other months of the year. All seasons have such periods. The QBO phase advances with these particular months and spans 21 to 29 such months in a cycle. The difference in the number of times a QBO phase occurs in a given month compared to another month is linked to the change of phase speed. The phase speed changes by month as well as by phase. When westerlies or easterlies blow at the upper vertical levels and winds of the opposite direction are present at lower levels, the phase speeds primarily fall between 0 and 0.1*π* a month through the year. In contrast, when easterlies dominate through the pressure levels, the phase speeds rise to 0.1 to 0.3*π* a month during April–November but fall between 0 and 0.1*π* a month during December–March. When westerlies prevail through the pressure levels, the phase speeds rise to 0.1 to 0.3*π* a month during September–June but mostly fall between 0 and 0.1*π* a month during July–August. A QBO phase is more recurrent in a month when its associated phase speeds tend to be zero and rarely occurs if its associated phase speeds are comparatively high. From the figures showing the QBO height–phase diagram (Fig. 9) and the monthly distribution of QBO phases (Fig. 11), the time of occurrence of phase transitions can be estimated. The figures also indicate that the easterly phase slows its descent speed during Northern Hemisphere winter.

At this point, it seems obvious to mention that supplementing the seven QBO winds with other variables (e.g., temperature field, meridional velocities, etc.) and conducting an NLPCA on this larger dataset will provide time series of nonlinear principal components that can be deciphered to investigate other physical connections within the QBO phenomenon.

## Acknowledgments

This study is supported by a grant from the Natural Sciences and Engineering Research Council of Canada and funding from the Canadian Climate Variability Research Network.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

### APPENDIX

#### Strategy of Model Optimization

The objective function (i.e., the mean square error of output to input signals) of the NLPCA model is minimized using a hybrid procedure (Webb and Lowe 1988) consisting of a quasi-Newton method for nonlinear optimization and a least squares method for linear optimization because the neural network model employs linear transfer functions for all the output neurons. The quasi-Newton method with a mixed quadratic and cubic line search procedure and the Broyden–Fletcher–Goldfarb–Shanno (BFGS) formula for updating approximations of the Hessian matrix (Broyden 1970a,b; Fletcher 1970; Goldfarb 1970; Shanno 1970) is used to minimize the objective function with respect to the weights and biases of bottleneck and decoding neurons. Every time these parameters are updated, the decoding signals are changed thereby, and a one-step exact minimization of the objective function with respect to the weights of output neurons is performed using the least squares method.

The quasi-Newton method is a fast algorithm to train feed-forward neural networks (Bishop 1995), but it stops at the first minimum it finds. The error function of nonlinear neural networks usually bears more than one minimum. As an effort to search for the lowest minimum, a model run is defined as follows: After a minimum of the objective function is found, the current values of the weights and biases will be randomly and slightly adjusted and then a new minimum is searched (Glover and Laguna 1998). If the new minimum is lower, the newly optimized parameter values will replace the current ones. Otherwise, the current values remain unchanged. The procedure repeats until no lower minimum has been found for several consecutive adjustments. Among all the runs of a particular neural network model, the solution with the lowest minimum of the objective function is taken as the model solution.

## Footnotes

*Corresponding author address:* Bei-Wei Lu, Department of Earth and Ocean Sciences, University of British Columbia, 6339 Stores Rd., Vancouver, BC V6T 1Z4, Canada. Email: blu@eos.ubc.ca