## Abstract

The concept of feedback has been used by several authors in the field of climate science to describe the behavior of models and to assess the importance of the different mechanisms at stake. Here, a simple 1D model of climate has been built to analyze the water vapor feedback. Beyond a static quantification of the water feedback, a more general formal definition of feedback gain based on the tangent linear system is introduced. This definition reintroduces the dynamical aspect of the system response to perturbation from Bode's original concept.

In the model here, it is found that, even though the water vapor static gain proves consistent with results from GCMs, it turns out to be negative for time scales below 4 yr and positive only for longer time scales. These results suggest two conclusions: (i) that the water vapor feedback may be fully active only in response to long-lived perturbations; and (ii) that the water vapor feedback could reduce the natural variability due to tropospheric temperature perturbations over short time scales, while enhancing it over longer time scales. This second conclusion would be consistent with studies investigating the influence of air–sea coupling on variability on different time scales.

## 1. Introduction

The concept of feedback gain was introduced by Bode (1945) to characterize the response of linear electric circuits to input signals. In spite of the nonlinearity of the climate system, the feedback concept is widely used in the climate science field to characterize the climate behavior, as in Coakley (1977), Cacuci and Hall (1984), Hansen et al. (1984), Wetherald and Manabe (1988), Cess et al. (1989), Andronova and Schlesinger (1991), Lindzen (1993), Schneider et al. (1999), Houghton et al. (2001), Bellon et al. (2003), and Friedlingstein et al. (2003), among numerous others.

To carry out such feedback analyses, models are privileged, because it is extremely difficult to extract such synthetic information on the climate behavior from observations, even though important details of the feedback processes can be gained from comparison between observations and simulations (e.g., Trenberth et al. 2005; Minschwaner and Dessler 2004; Allan et al. 2003; Soden et al. 2002). However, the practical determination of feedback gains from complex simulation models is far from direct and, as we emphasize in this paper, the feedback concept is meaningful only if one can properly define the feedback loop in the model.

Hansen et al. (1984), among others, analyzed climate sensitivity in terms of leading feedback mechanisms. They use the following definition: when a feedback loop going through the variable *ϕ* has been defined, a feedback is characterized by its gain (*g*) or its response factor ( *f* ), which is defined by

where *δϕ*^{∞} is the change in equilibrium value of *ϕ*, after a forcing perturbation has been applied; *δϕ*^{0} is the change in the equilibrium value of *ϕ* for the same perturbation but when the feedback is suppressed.^{1} This usual feedback gain definition is, therefore, based on the difference between two equilibrium values, and will henceforth be called static gain.

Another approach, used by Coakley (1977) and Wetherald and Manabe (1988), has the advantage of detailing the intervening elementary processes within the feedback loop. To assess the role of a set of feedbacks in the response to a change (Δ*F*) in the net radiation at the top of the atmosphere (*F*), they use partial derivatives of the outgoing flux (*R*) with respect to the variables involved in each feedback (*x _{i}*,

*i*= 1, . . . ,

*n*) and the full derivatives of each involved internal variable with respect to the mean surface temperature (

*T*). Because incoming and outgoing fluxes are equal at equilibrium, a relationship between the incoming flux perturbation and the equilibrium mean surface temperature is obtained:

_{s}One can notice that by dividing this equation by ∂*R*/∂*T _{s}*, Eq. (2) can be rewritten in the form of Eq. (1):

This method has also been extended by Colman et al. (1997) to account for second-order derivatives in order to give an idea of the nonlinearity of the processes.

As a first constraint on practical applications of this formulation, the perturbation of the system has to apply primarily to *T _{s}*, and must not directly influence any other variable. Hence, in the case of climate change, considering a CO

_{2}concentration perturbation requires the assumption that a CO

_{2}change impacts only the surface temperature (which may modify the lapse rate) but does not directly impact the lapse rate. A second implication of this definition concerns the summing of individual feedback gains to determine the complete response of climate. The full derivative factor

*dx*/

_{i}*dT*in Eq. (3), which determines the individual gain

_{s}*g*, must not include processes that are in common between different

_{i}*g*terms. Practically, however, it may be difficult to separate these different components in the model response.

_{i}Additionally, applying a linear concept to nonlinear climate models requires some assumptions: it is necessary that, as soon as all nonlinear rapid transient processes have been damped out, the most inert processes of climate respond linearly to perturbation of moderate amplitude. Recent results obtained by Goodman and Marshall (2002) suggest that this assumption might be correct.

The various difficulties raised by feedback definition can be summarized in three problems: (i) a structure problem, because each feedback loop must be defined unambiguously in the system, (ii) the choice of a perturbation adequate to a proper definition of gain (e.g., a SST perturbation or a tropospheric temperature perturbation), (iii) the nonlinearity of the model.

Another difficulty in the interpretation of feedback processes arises from the speed of the different responses. Some processes participating in the feedback mechanism may be fast, others may be very slow. This problem is usually avoided by analyzing only the new equilibrium reached by the system after a perturbation has been applied. This current practice might, however, ignore important dynamical components of the response, especially when the forcing is varying.

To investigate these issues and illustrate our methodological approach, a simple 1D model for the atmospheric water vapor (WV) is introduced within a particular formalism (section 2). We show in section 3 how the usual feedback concept can be extended in order to include the dynamical characteristics of the climate response. In section 4 the model is used to characterize the dynamics of the WV feedback, and we derive some results on the atmosphere dynamics and variability. The results are summarized in the final section and remaining difficulties concerning the definition of feedback loop in complex systems are discussed.

## 2. Model

The model is intended to focus on the global WV dynamics, with physical processes other than radiative transfer being crudely reproduced. After having presented the main assumptions, we introduce our specific formalism, which will be used to generalize the definition of feedback.

### a. Description

A single column of atmosphere with three cloud layers is considered. Figure 1 displays a schematic diagram of our model. Crude assumptions are applied: (i) the fact that we use a single-column model rules out the possibility of any dynamical feedback, such as changes in the Hadley circulation, midlatitude eddies, etc., that would affect water vapor feedback; (ii) convection is not explicitly modeled but its effects are taken into account through fixed lapse rates. This can be justified by Zhang et al. (1994) and also Colman (2003), who showed that variation in lapse rate does not alter significantly the combined WV and temperature lapse rate feedback; (iii) ocean mixed-layer depth is fixed; (iv) stratosphere shortwave absorption is fixed; (v) no explicit cloud-cover modeling is introduced, for simplicity reasons and because cloud feedback is not robust even in GCMs (see, e.g., Held and Soden 2000; Colman 2003).

Atmospheric water vapor content is driven by evaporation and precipitation. Evaporation (*E*) depends classically on the saturation specific humidity at the ocean mixed-layer temperature and on the surface relative humidity (RH), the drag coefficient being fixed. The surface RH *r*^{s}_{h} is calculated using the column water vapor content *Q*, the troposphere temperature *T _{τ}* and a fixed relative gradient of RH, which is such that RH at tropopause is half the surface RH:

where *H*_{ttp} is tropopause altitude.^{2}

The column water content and the surface RH are then linked by the relationship

where *P*_{sat}(*T*) is the water vapor pressure at saturation (Clausius–Clapeyron), *H _{m}* is the altitude of mean troposphere,

*T*is the tropospheric temperature at altitude

_{τ}*H*, and

_{m}*γ*is the tropospheric lapse rate.

_{τ}Deep convection is the main mechanism constraining RH in tropical regions. In extratropical regions, RH is driven by synoptic perturbations. We represent very roughly these processes by modeling the precipitation amount (*P*) as equal to evaporation plus a term driving RH toward a target profile of RH, *r*^{∞}_{h}(*z*), which is supposed to be fixed independently of the CO_{2} concentration. This target profile of RH is defined by a target surface RH, *r*^{s,∞}_{h} = 0.7, and Eq. (4). Then, we can calculate a target water vapor content *Q*^{∞}(*T _{τ}*) using Eq. (5). Precipitation

*P*is then given by

The symbol *τ _{p}* refers to the water content return-to-equilibrium characteristic time. Betts and Miller (1986) use a similar humidity target for precipitation modeling, and our model is chosen to be consistent with results on longer time scale from Hall and Manabe (2000a), who showed that precipitation annual mean is driven by evaporation annual mean, and with the observation that RH is near constant at equilibrium in GCMs, whatever the CO

_{2}concentration change—see Houghton et al. (2001), chapter 7, Salanthé and Hartmann (1997), Del Genio et al. (1994), or Hansen et al. (1984)—or the shortwave (SW) incoming flux change (e.g., Soden et al. 2002).

A strong assumption of this model is that precipitation is driving RH toward a target value, whatever evaporation is. Hence, when an evaporation perturbation occurs, it leads to an instantaneous change in precipitation, without any change in atmosphere water content. Also, an increase in troposphere temperature, at fixed evaporation, leads to a decrease in precipitation because of the corresponding stabilization of the atmosphere. This is consistent with the fact that the effect of an increase in CO_{2} concentration without surface latent flux adjustment is a decrease of precipitation (e.g., Yang et al. 2003). Precipitation then adjusts back to evaporation over the time scale *τ _{p}*.

The target RH is set to 0.7 at surface level, and the characteristic time of return to equilibrium is set to 30 days, as a compromise between the synoptic characteristic time (a few days) and the radiative time scales, or Hadley cell overturning times, that can reach up to 30–70 days. The sensitivity of our results to this important model parameter will be analyzed in section 4.

The radiative module is a 65-layer column model of the atmosphere, with three cloud layers and two gases (H_{2}O and CO_{2}). The module computes longwave radiative budgets of troposphere, stratosphere, and ocean, using a Malkmus narrow-band model with a water vapor continuum. The principles behind this module have been explored by Green (1967) and developed later by Cherkaoui et al. (1996). Radiative modules based on these principles were built by Hartmann et al. (1984) and Soufiani et al. (1985), who also created the radiation band coefficient tables.

To allow for dynamic feedback analysis, we need to introduce the Tangent Linear System (TLS) of the model. For this purpose, we make use of the transfer evolution formalism (TEF, described below) to mathematically model our system.

### b. The transfer evolution formalism prescriptions

The TEF is a mathematical framework for system analysis and simulation [see appendix A or our Web site (http://www.lmd.jussieu.fr/zoom) for a more detailed description]. TEF is appropriate to our problem since it systematically uses the TLS of a model and, as will be shown later, the TLS is central to the feedback gain definition. The model presented in previous section is mathematically represented by a set of equations corresponding to two kinds of objects:

- cells, which are elementary models and correspond to state equations such as
is the state variable vector of cell**η**_{α}*α*and {*ϕ*} represents dependent boundary conditions; that is, variables considered as boundary conditions by a cell, depending upon the complete model state. These dependent boundary conditions are required to make each cell correspond to a well-posed problem. (The five model-state equations are described in Table B1.)_{i} - Transfers (described in Tables B2 and B3), which are determined by static constraint equations such as As a simplified notation, let
also be the state vector of the complete system and*η*be the vector of the full dependent boundary conditions. When initial conditions are given at time*ϕ**t*_{0}, the system is a well-posed problem.

To solve the system at each time step, the TLS of the model is derived around its current state^{3} [** η**(

*t*)] by writing

_{n}**(**

*η**t*) =

**(**

*η**t*) +

_{n}**(**

*δη**t*) and

**(**

*ϕ**t*) =

**(**

*ϕ**t*) +

_{n}**(**

*δϕ**t*), giving

This method of obtaining the TLS has been formally introduced by Cacuci (1981a, b. We use the notation *δ̃* to make explicit that the system (9) is a second-order approximation in *δt* = (*t* − *t _{n}*) of the nonlinear system defined by Eqs. (7) and (8).

Such a linear system can be solved via a symbolic transformation method. We prove in appendix A that the Borel transform^{4} of this TLS can be written as

where ℬ[ *f* ](*τ*) is the Borel transform of *f* (*t*), *τ* is the Borel variable, ** δ̃η**(

*τ*) and

**(**

*δ̃ϕ**t*) are the solutions of the TLS, and where the quantities

*δ̃η*_{dec},

**F**,

**C**,

*δ̃ϕ*_{ins}can be calculated from elementary Jacobian matrices and vectors at time

*t*.

_{n}All numerical results are obtained with software developed by the present authors and colleagues to implement models expressed with the TEF (more information is available on http://www.lmd.jussieu.fr/zoom). An approximation of the step-by-step evolution of the complete system is obtained by solving system (10) thanks to the approximations

for small *δt*.

### c. Validation of the model

The TEF implementation of the model results in five cells and seven transfer objects, with a total of 18 variables. Some basic tests are now reported to assess the realism of this simple model concerning the WV feedback. The change in equilibrium is first considered, then the fast and slow responses to CO_{2} concentration changes are discussed. At equilibrium, with a CO_{2} concentration of 300 ppm, the model gives the following results. The SW fluxes absorbed by the different objects are 40 W m^{−2} by stratosphere, 16 W m^{−2} by troposphere, and 164 W m^{−2} by ocean. In comparison, Salby (1996) gives 60 W m^{−2} for the atmosphere (troposphere plus stratosphere), and 170 W m^{−2} for the ocean. The longwave (LW) flux exchanged in the model by the different components of the climate system at the 300 ppm equilibrium are also close enough to the estimates from observation (e.g., Salby 1996, p. 43) to assume that the LW radiative module is able to reproduce realistic budgets: the direct emission from surface to space is about 15 W m^{−2}, to be compared with 20 W m^{−2} according to measurements; the exchanges between surface and atmosphere give a positive flux of 41.1 W m^{−2}, to be compared with 50 W m^{−2}; and the outgoing longwave radiation (OLR) is 221 W m^{−2}, to be compared with 220 W m^{−2}. It should be mentioned that the other exchanged fluxes (e.g., sensible flux at the surface of 24 W m^{−2} or latent heat flux of 83 W m^{−2}) are also close to the values found with complex models, showing that the basic equilibrium budgets are matched within the simple model.

Figure 2 shows the time evolution of the state variables, in response to a step from 300 to 600 ppm in CO_{2} concentration. On the stratospheric time scale, troposphere and ocean temperatures are almost constant. Consequently, the stratosphere reaches a colder quasi equilibrium corresponding to the new CO_{2} concentration and to the tropospheric and oceanic initial temperatures. This situation corresponds to the standard definition of CO_{2} radiative forcing. The results from the model (e.g., at the doubled CO_{2} concentration, *F*_{2}* _{X}* = 3.3 W m

^{−2}) are fairly close to those of other studies (Houghton et al. 2001, chapter 6).

On longer time scales, troposphere and ocean warm up slightly, and stratosphere, always at radiative quasi equilibrium, begins to warm up slowly. Nevertheless, the final effect on stratosphere is a strong cooling, since the initial stratosphere cooling is much larger than the subsequent warming.

The climate sensitivity of the model to doubling CO_{2} concentration is +1.3 K, with an increase in water vapor content of 3.3 km m^{−2}, corresponding to a 10% rise. This is low compared with most GCM results (see Houghton et al. 2001, chapter 9; Kothavala et al. 1999), but is still a classical value for 1D models.

## 3. Methodology for feedback study

### a. Definition of the water vapor feedback

Within the TEF, a feedback loop is defined as a set of processes interfaced by transfer (vector) variables {*ϕ _{j}*,

*j*= 1, . . . ,

*n*} in which the evolution of each variable

*δϕ*depends only on

_{j}*δϕ*

_{j−1}(and the evolution of

*δϕ*

_{1}depends only on

*δϕ*). This definition is consistent with the terminology recalled in section 1 for static gain.

_{n}Our simple model is expressed in order to implement one among the possible definitions of WV feedback. The loop in the model is the following: a troposphere temperature increase occurs; RH is decreased (as a consequence of the Clausius–Clapeyron relation); precipitation, modeled as to maintain a target RH, decreases; RH goes back to its initial level, corresponding to a larger total water content; the radiative budget is modified (the troposphere and ocean warm up); as a consequence, the troposphere temperature increase is amplified. The system is supposed to be perturbed by a change in CO_{2} concentration: the perturbation hence applies to the tropospheric temperature first, because of the modified radiation budget.

The first obstacle to analyze the WV feedback is the fact that, in our model, a WV feedback loop cannot be defined as a set of variables and processes verifying the required conditions. For instance, troposphere temperature (*T*) influences many processes that do not belong to the WV feedback (e.g., the LW fluxes depend on *T* independently from the WV feedback). In other words, *T* is involved in several feedback loops so that the feedback gain associated with the *T* variable results from an interplay between many feedbacks.^{5} To define a feedback loop, we need at least one variable that is involved only in this loop.

To obtain loops that involve mechanisms that can be separated, one might have to explicitly reshape a model with an unambiguous structure. In the present case, the model is modified by introducing an extra variable: *T*_{WV}, defined as the temperature driving the water-vapor-related mechanisms. In a more complex model, *T*_{WV} might be the temperature at a given atmospheric level or a temperature modified by some processes. In the present model, *T*_{WV} is taken as equal to the atmospheric temperature *T*. In other words, the unperturbed *T*_{WV} is set equal to *T* but does not represent the same physical parameter: a perturbation can apply to *T*_{WV} independently of *T*. The introduction of this variable modifies the model equations in the following way:

where *T*, *Q* are the mean troposphere temperature and humidity, respectively, and *η* represents all the other variables. The functions *f*, *F,* and *G* correspond to the model equations. Although the two models have identical solutions for *T*, *Q*, and *η*, the feedback loops are different (see Fig. 3). In the new structure of the model, there is a feedback loop going through *T*_{WV} that is distinct from all other loops.

This methodology, akin to the one used by Hall and Manabe (2000b) in their GCM feedback analysis, enables one to suppress the WV feedback without losing water conservation and model consistency.

### b. Methodology for dynamic characterization of the water vapor feedback

The TEF algorithms solve the TLS to provide the full dynamics of the perturbed system. When the system is nonlinear, the TLS may evolve with time. We will illustrate the concept using the asymptotic stable equilibrium state of the system where the TLS is autonomous. Doing so, we leave aside the important problem of dealing with nonlinear systems out of equilibrium.

Once the TLS of a model is available in the form of Eq. (10), the method to characterize a specific feedback gain in a system is straightforward: it consists in eliminating all variables except one—here *δ̃T*_{WV} for the WV feedback—in the algebraic system (10). The remaining scalar equation reads

It is shown in appendix A that *δ̃T*_{WV,ins}(*t*) is the variation of *T*_{WV} when the feedback loop is cut just after *T*_{WV} in Fig. 3. The function *δ̃T*_{WV}(*t*) is the variation of *T*_{WV} in the complete TLS.

We define *g*_{WV}(*τ*) as the dynamic gain of the WV feedback because it links the dynamics of the full model [*δ̃T*_{WV}(*t*)] to the dynamics of the open loop model [*δ̃T*_{WV,ins}(*t*)]. It thus represents the effect of closing the feedback loop: Δ*T*_{WV} → *impacts on the rest of the system* → *feedback perturbation on T*_{WV}. It is noteworthy that the calculation of the feedback gain is here realized without actually applying any perturbation to the system, but instead is carried out semianalytically from the very model equations.

Note that function *g*_{WV}(*τ*) generalizes the concept of static gain, since at the infinite limit

While the static feedback gain gives only the response corresponding to the new equilibrium value, the dynamic gain describes in addition the whole response dynamics of *δ̃T*_{WV}(*t*), and hence the whole dynamics of the feedback, through the inverse Borel transform of Eq. (13) rewritten as

where the “*” sign stands for the convolution product.

The time-varying function ℬ^{−1}{[1 − *g*_{WV}(*τ*)]^{−1}}(*t*) × Δ*T*_{0} may be interpreted as the change in *δ̃T*_{WV} after a Δ*T*_{0} step is applied on the *T*_{WV} equation, corresponding to the following modification in the equations of (12): *T*_{WV} = *T* if *t* ≤ 0, and *T _{WV}* =

*T*+ Δ

*T*

_{0}if

*t*> 0. This response includes the perturbing step, so that

*δ̃T*

_{WV}(

*t*) is not continuous at

*t*= 0.

The response function can now be defined as

From Eq. (15), and remembering that the Dirac distribution *δ*(*t*) is the derivative of the unit step function ϒ(*t*), *F _{TWV}*(

*t*) has two different interpretations: (i)

*F*(

_{TWV}*t*) is the response of

*T*

_{WV}within the full model, when the

*T*

_{WV}equation is perturbed by a unit step

^{6}ϒ(

*t*); (ii)

*F*(

_{TWV}*t*) is the response of

*T*

_{WV}when feedback is active, to a perturbation that would have led to a unit step response if the feedback was suppressed.

### c. Characteristics of the water vapor feedback in the model

Thanks to explicit calculations of all Jacobian matrices in our model, *g*_{WV}(*τ*) can be computed for a range of values of the Borel variable (*τ*). Then, the Borel transform of the feedback function can be identified using a fit to elementary transform functions, as explained in appendix A. Applying this method, we found that two elementary functions fit the response function in the Borel space to an accuracy of 10^{−6}. The corresponding function in the physical space could be determined as

with *τ*_{1} = 36 days and *τ*_{2} = 4.8 yr, when *τ _{p}* = 30 days is chosen in Eq. (6).

Because the value of *τ _{p}* is uncertain, a sensitivity analysis is carried out where

*τ*is varied from 5 days to 1 yr. The wide range on

_{p}*τ*is given for completeness, but its value should not exceed a few months. As an indication of this, Hall and Cacuci (1983) give a value of 27 days as the time of return to equilibrium of the perturbed atmospheric temperature along a trajectory that includes convective events. The results are reproduced in Table 1 and the corresponding response functions are shown in Fig. 4.

_{p}The first pole characteristic time is close to the precipitation characteristic time *τ _{p}*, but the second is much less sensitive to it. When

*τ*is longer, that is, when RH is less nudged to its equilibrium value, the negative part of the WV feedback becomes less important. This negative part, however, is still visible for the largest values of

_{p}*τ*. One can also notice that the two residues are quite insensitive to

_{p}*τ*.

_{p}## 4. Discussion of the results

Our simplified model has been used to illustrate how the TLS can base a rigorous definition of feedback functions and allow for a precise analysis of the dynamics of the perturbed system. The model was also built to get a better understanding of the dynamics involved in WV feedback processes, and we now discuss the realism of the results. The origin of the two poles is first analyzed, and the response to realistic perturbations is given. An application to climate variability serves also as a comparison with the results of Barsugli and Battisti (1998).

### a. Interpretation of the feedback function

The fast pole corresponds to the lowering of latent heat flux due to rainfall decrease, which comes from the rising temperature (corresponding to a decrease in RH). This mechanism constitutes one path of the WV feedback: any transient trajectory showing an increase in atmospheric absolute humidity requires an imbalance between precipitation and evaporation, and hence necessitates an increase in atmospheric latent energy content compared with the equilibrium state. In consequence, the WV feedback process should involve a rapid atmospheric cooling, as formalized in our model, with a time response of about a few weeks. This negative feedback can be illustrated by considering a doubling CO_{2} experiment that would lead to a 3 K surface temperature increase: for constant RH, absolute humidity would increase by about 10%. This would correspond roughly to an addition of 4 kg m^{−2} of water vapor in the atmosphere, and hence to a latent energy loss of 10^{7} J m^{−2}. Because the involved LW flux changes are about 3 W m^{−2}, one month of the total additional flux would be necessary to collect that energy.

Barsugli and Battisti (1998, hereafter B98) analyzed the effect of atmosphere–ocean coupling on midlatitude variability using a model similar to ours: no cloud feedback, constant RH, ocean mixed layer of 50 m depth, and perturbation applied to atmospheric temperature. The main difference with our model resides in their gray gas radiative model. Their model is a constant coefficient TLS of two equations for the troposphere and ocean temperatures anomalies. The model is analyzed in the Fourier space, and analytic solutions are computed. We applied the Borel transformation to their model and found a slow pole of about 6 months with positive feedback, and a fast pole of about 4 days with negative feedback. The amplitude of the response associated with their negative pole is −0.72 in B98, whereas Table 1 gives −0.35 in our model. Their results are thus consistent with ours, even though not identical.

The fast feedback mechanism can also be associated with the results of Frankignoul et al. (2004) for whom the short-term negative feedback found in GCM experiments is due to atmosphere–ocean coupling through latent and sensible heat fluxes. Unfortunately, they did not discuss the influence of precipitation on the feedback value.

The slow pole corresponds to the classical WV feedback, as being the process by which atmosphere warms up when an increase in water vapor concentration modifies the LW radiative balance. Table 1 gives an estimate of 4 to 7 yr as characteristic time of that slow pole.^{7} This very long building up time of the feedback can be explained by the fact that WV feedback is mediated by successive nonzero characteristic time processes: an initial tropospheric temperature perturbation leads to an increase in absolute humidity; then, the consecutive change in the radiative fluxes increases progressively the SST, leading to a new change of the tropospheric temperature, which will, in turn, reenter the feedback loop. This step-by-step process leads to a final temperature increase, involving long characteristic times and complex behavior, ignored by the equilibrium approach of static feedback analyzes.

Our characteristic time of 3 to 7 yr is very long compared to the 6 months found with B98. Considering the very crude radiative model used in B98 (a gray gas model), we take this as a confirmation that the long response time of the positive feedback gain is linked with the radiative interactions between ocean and atmosphere.

These different features of the response function show that (i) it is essential to consider the whole feedback loop to assess the characteristic time of a feedback process: even if RH responds very quickly to temperature changes, the WV feedback is very slow; (ii) the WV feedback in this model is fully active only in response to a perturbation that lasts more than ten years, making its assessment with short-term experiments (like volcano eruptions) particularly difficult, as we will discuss later.

The WV static gain can be retrieved by applying Eq. (14) to Eq. (17), yielding a value of *g*^{stat}_{WV} = 36%, which is close to the results from GCMs (see Lindzen 1993; Schneider et al. 1999). It is also close to 0.40, the value given from theoretical considerations by Held and Soden (2000) in their review article.

Nevertheless, the elicitation of the fast and negative part of the feedback demonstrates the interest of the methodology: the WV feedback is not a monotonous process only responsible for enhancing the temperature increase but rather, a complex dynamic process responsible for a temporal pathway of the temperature change. Some consequences of that understanding of WV feedback are discussed now.

### b. Response of the model to a realistic perturbation function

The choice of the unit step function is made to ease the mathematical definition of the feedback function. Equation (15) provides, however, the response of the model to any perturbation. As an example, Fig. 5 represents the response of the model with WV feedback to a perturbation that would lead in the no-feedback model to a ramp response^{8} of 0.1 K decade^{−1}, beginning at *t* = 0. It can be seen that the negative feedback effect is hardly detectable but is responsible for a reduction of the warming during the first 5 yr. The total response starts with the influence of the fast pole only, and its slope is progressively built up with the influence of the slow pole.

### c. Influence of a feedback pole on variability

Since the feedback is now dynamically characterized, it is possible to endow the model dynamics with the influence of each feedback pole. This can be made explicit by considering a sinusoidal perturbation with amplitude *A*_{0} and pulsation *ω* starting at *t* = 0, such that the temperature variation is given by Eq. (18) if the feedback is made inoperative:

As shown by Eq. (15), the overall temperature response to this perturbation is given by the convolution product:

Assuming that all poles are real, negative, and simple, the response function *F _{TWV}*(

*t*) reads

Elementary calculations lead to

where 1) is the perturbation, 2) is the static influence of each pole,^{9} 3) is the transient influence of the pole. This equation shows that, for each feedback pole *i:*

for high frequencies, when

*ω*^{−1}≪*τ*, the pole has no influence on temperature variation;_{i}- for low frequencies,
*ω*^{−1}≫*τ*, the pole does have an influence on temperature variation. In that case, the influence of the_{i}*i*th pole leads to the following permanent temperature response: Hence, if such a feedback pole corresponds to a positive feedback, its effect is to enhance the oscillation amplitude. Conversely, a negative feedback pole reduces the variability. Between these extremes, when the feedback pole characteristic time is close to the perturbation time scale, the influence of the pole decreases as the perturbation frequency increases.

Taking now *A*_{0} = 1 and dropping the third rhs transient term, Eq. (21) gives the variance ratio between the responses of the system with and without WV feedback to tropospheric temperature perturbations. Assuming that natural variability comes, at least partly, from such perturbations, we can use this equation to assess how these perturbations will be amplified or damped by the feedback.

For high frequencies, only the fast negative pole of the WV feedback is active, and the feedback is found to reduce natural variability. For low frequencies, on the other hand, both poles are active, and the WV feedback is found to eventually increase natural variability. These results illustrate how a time-varying feedback can influence the power spectrum of natural variability.

Our variance ratio can be compared with the power spectra of atmosphere temperature in the B98 model with and without air–sea coupling (B98; Fig. 4), that shows that air–sea coupling increases variability at low frequencies and reduces it at high frequencies. In our model, the variance ratio is becoming smaller that unity for frequencies greater than about 1.4 × 10^{−4} cycles per day (cpd), when B98 found this frequency to be approximately 7 × 10^{−3} cpd. This value is pulled toward slower frequencies in our case because of the influence of the slow pole: we checked that when our slow pole is set to 6 months, we retrieved B98 results.

### d. Evidence of the fast negative feedback gain in literature

In our model, the fast negative feedback pole is primarily due to a negative coupling between precipitation and air–temperature perturbation, arising from Eq. (6) and the Clausius–Clapeyron relation. This model feature is also observed in GCMs, as shown in Yang et al. (2003): in their study, the authors forced a GCM by an instantaneous doubling of the CO_{2} concentration, and observed a reduction of convection, evaporation, and precipitation over short time scales (a few days, see their Fig. 1). This reduction of the hydrological cycle transfers a large part of the additional energy from atmosphere to ocean, increasing the SST and cooling the atmospheric temperature. Over long time scales (a few years), a new warmer equilibrium is reached, in which a higher SST overcomes the reduction of convection, leading to a more intense hydrological cycle only when the new equilibrium is reached.

The instantaneous CO_{2} doubling experiments are informative on the climate system, but are unrealistic. As for the real world, it is extremely difficult to find in published observations a confirmation of the existence of short term anticorrelation effects, mainly because correlation studies use SST anomalies for practical reasons. According to our results, the anticorrelation should be only observed when these SST anomalies are driven by air–temperature anomalies. Kang et al. (2004) observed such a negative correlation in the subtropical western Pacific, a region where air–sea interactions are active. In opposition, their results from a seasonal prediction system forced by the SST (i.e., without air–sea feedback) shows none or even slightly positive correlation. They concluded that short-term air–sea coupling was at the origin of the anticorrelation, and we have shown that this could be related to the existence of the negative feedback pole.

Several model-based analyses found the air–sea coupling influence on variability dependent on time scale. Barsugli and Battisti (1998) is one of those, as mentioned in the previous section. In a subsequent paper, Bretherton and Battisti (2000) compared two GCMs with a stochastic model derived from B98. They demonstrated that the conclusions drawn from B98 were consistent with analyzes of the North Atlantic Oscillation in GCMs. Another study from Bhatt et al. (1998), using the National Center for Atmospheric Research (NCAR) Community Climate Model (CCM1), reached the conclusion that “air–sea interaction has a greater impact on seasonally averaged variance than monthly variance”, which again is consistent with the existence of a short-term negative feedback and its consequences on variability, as described in section 4c.

The climate response to the Mount Pinatubo eruption, which has been comprehensively measured and analyzed (e.g., Robock 2000, 2002), has been used by several authors to validate models and/or to assess the climate sensitivity (e.g., Soden et al. 2002; Forster and Collins 2004; Wigley et al. 2005). However, during volcano eruptions, the perturbation is mainly applied to the SW flux absorbed by the ocean. The ocean, and its large inertia, is therefore driving the change, and atmosphere is slaved to the ocean. In such a situation, the atmospheric temperature increase is a ramp and, as can be seen in section 4b, the negative pole of the feedback is only responsible for a slight reduction of the slope of the tropospheric temperature change, explaining that its detection is much more difficult.

Concerning the time scale of the whole feedback, despite the fact that Soden et al. (2002) show clearly that the water vapor feedback is active during the few years following the Pinatubo eruption, it does not follow that this feedback is fully active in our sense: according to Fig. 4, *F _{TWV}* lies between 1.2 and 1.4 at

*t*+ 5 yr; that is, the water vapor feedback is found to increase the warming (or, in this case, the cooling) by between 20% and 40% after 5 yr, to be compared with 60% at the infinite limit.

## 5. Concluding discussion

The analysis presented here uses the tangent linear system of models as the fundamentals of the feedback gain definition. Following this method, we built a simple model including a detailed radiative transfer scheme to analyze water vapor feedback processes. We have shown that even if the static gain is unchanged, our definition introduces the dynamic features that recover Bode's original definition based on the Laplace transformation. The dynamic feedback functions of the Borel variable *τ* unveil important characteristic times of climate mechanisms.

In our model, the WV feedback is found to have a positive static gain of 36%, and a characteristic time longer than 4 yr, making the WV feedback fully active only in response to perturbations that last at least 10 yr. Except for the longer response time due to our more detailed treatment of the atmospheric humidity and radiative exchanges, our model give results consistent with those of Barsugli and Battisti (1998).

The more remarkable result is that our assumption regarding relaxation to a fixed RH target can lead to a transient negative WV feedback, which extends for far longer than the associated RH relaxation time. These dynamic features are well recognized because the perturbation needed by the mathematical definition of feedback gain is a unit step function. For realistic perturbations, the responses are determined by a convolution product that may fade out the short-term component of the feedback function, making them difficult to observe.

We have shown, however, that the poles of the response function are taking part in the determination of the variance ratio between coupled and uncoupled systems, as studied by Barsugli and Battisti (1998). Particularly, our findings suggest that WV feedback might reduce short-term natural variability in response to troposphere temperature perturbations, and, on the other hand, increase the long-term natural variability.

This analysis emphasizes that complex interactions between processes—even if it is limited here since it concerns a single-column model—lead to a dynamical complexity of the climate response. We have also shown that it is far from trivial to separate feedback loops in a strict manner, even in simple models. This separation is essential to justify the classical methods, which retrieve feedback gains from GCMs, a difficulty already reported by Schneider et al. (1999). We proposed a methodology to do so, which is based on the introduction of extra variables to separate the various mechanisms.

Although our model is nonlinear, the analysis was entirely based on the TLS at equilibrium. There is no obvious way to extend the method to nonlinear systems out of equilibrium, for which the TLS becomes nonautonomous. Intuitively, the long-term response function should not be sensitive to the variation of the TLS along the system trajectory, as long as the numerous processes that are not represented in our simple model (e.g., large-scale circulation, thermohaline circulation, sea ice cover, . . .) are not too much impacted. This thinking meets the conjectural conclusions of Goodman and Marshall (2002) who compared the singular elements of the TLS of their quasigeostrophic model with the EOFs of the nonlinear long term simulations. They found close analogy between the slowest singular modes and the leading EOFs, suggesting that the slowest structures are responding linearly to the quickest (nonlinearly filtered) synoptic perturbations. As a consequence, the linear slowest response should reproduce a realistic long-term climate behavior, as long as no bifurcation occurs. If the slowest modes are robust, even when the fastest ones are uncertain, then the inert part of TLS might be quasi autonomous.

An original direction of research is given by Aires and Rossow (2003), who characterize nonlinear systems following the nonautonomous TLS as it evolves along a full model trajectory. In that case, the nonlinear system is described as histograms of the TLS transition matrix coefficients sampled along the trajectory. It is theoretically possible to extend this method with a dynamical characterization of each TLS coefficient in the histograms, but this leads to severe difficulties and there is still a long way to go before a comprehensive treatment of the climate system can be achieved.

Our own direction of research retains a combination of both previous methods. One can recover the feedback functions behavior *f* (*T*, *t*) along the trajectory between a perturbation time *t* and a probe time *T*, in analogy with the Borel variable *τ* = *T* − *t*, and verify that for *t* ≪ *T*, the functions become insensitive to *t*, linking the slow response to perturbation with Lyapunov exponents.

## Acknowledgments

The authors wish to thank Olivier Tallagrand for helpful comments and remarks, and Robert Franchisseur for computer support. Stéphane Hallegatte was funded by the Centre National de Recherches Météorologiques, Météo-France, and by the Ecole Nationale des Ponts-et-Chaussées. Three anonymous reviewers are thanked for their useful comments and suggestions.

## REFERENCES

**,**

**,**

**,**

**,**

**.**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

_{2}O.

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

_{2}and H

_{2}O applied to radiative properties and conductive-radiative transfer.

**,**

**,**

**,**

**,**

**,**

**,**

### APPENDIX A

#### The Transfer Evolution Formalism

##### Tangent linear system analysis

As explained in the article, the model is mathematically represented by a set of equations of two kinds:

Let ** η** be the state vector of the complete system and vector

**be the dependent boundary conditions. With initial conditions in**

*ϕ**η*at time

*t*

_{0}, the system is a well-posed problem.

The method consists in building the first-order development of the dynamical system around its current state [*η*(*t _{n}*)]. For each cell

*α*, it reads

where * δη_{α}*(

*t*) =

*(*

**η**_{α}*t*), and

_{n}**(**

*δϕ**t*) =

**(**

*ϕ**t*) −

**(**

*ϕ**t*).

_{n}The tangent linear system corresponding to system (A3) is, for each cell *α*,

where the suffix *β* sweeps the list of cells.

We approximate the true time evolution of the model [* δη_{α}*(

*t*) and

**(**

*δϕ**t*)] by

*δ̃η**(*

_{α}*t*) and

**(**

*δ̃ϕ**t*), the TLS solutions, since they differ only by 𝒪[(

*t*−

*t*)

_{n}^{2}].

In formulation (A4), the Jacobian matrices contain critical information for the analysis of interaction between variables. The TLS can be solved by various methods, including Laplace transforms. Rather than Laplace, we shall use the more convenient Borel transformation defined by

where *f̃*(*p*) stands for the Laplace transform of *f* (*t*). Contrary to the Laplace variable, the Borel variable *τ* is real and homogeneous to time.

Because ℬ[∂*f* /∂*t*](*τ*) = (1/*τ*)ℬ[ *f* ](*τ*), the Borel transform of Eq. (A4) reads

If the cell variables ** δ̃η** are eliminated from the second equation, the complete system of equations (which includes cells) becomes

where the quantities ℬ[*δ̃η*_{dec}], **F**, **C**, ℬ[*δ̃ϕ*_{ins}] depend on *τ* and can be calculated from the elementary Jacobian matrices and vectors at time *t _{n}*.

The first equation of (A7) describes the evolution of state variables. State variables evolve because (i) of their internal inertial evolutions *δ̃η*_{dec} (which would be obtained if transfer models were changed to constant model with ** δ̃ϕ** =

**0**), and (ii) of the evolution of their boundary conditions (

**≠**

*δ̃ϕ***0**). Matrix

**F**describes the influence of transfer variables on state, and independently of the type of model used for these transfers (

**F**is independent of the model of

**).**

*δ̃ϕ*In Eq. (A7), *δ̃ϕ*_{ins} represents the variation of transfer variables if ** δ̃η** =

*δ̃η*_{dec}(i.e., if the cell models were changed to decoupled models with

**F**= 0). Consequently,

**C**represents the effect of cell and transfer coupling.

The developed expression of matrix **C** shows how the partial derivatives defined at the cell and transfer level combine. The coefficients of the coupling matrix are rational fractions of the variable *τ*. This is the way the full dynamic of the system bounds the remaining variables after an algebraic elimination process.

##### Numerical solution of the transfer evolution formalism

For large systems, the above matrices are huge and sparse, and exhibit an internal structure that depends upon the connections between cells and transfers. The full algorithm of the ZOOM solver (ZOOM is a TEF-dedicated solver developed by the authors and colleagues) follows a technique called “relaxed super-nodes hyper multi-frontal method” (cf. Liu 1992). We focus here on the principles of the resolution that explain how the system dynamics is described by the coupling coefficients.

###### Equivalence between Borel transform and the Crank–Nicolson scheme

It is easily shown that the Crank–Nicolson resolution of system (A4) with a time step *δt*, is identical to its Borel transform (A7), with the correspondence *τ* ↔ (*δt*/2). To demonstrate this equivalence, let *δ̂X* be the time evolution of variable *X* approximated by a Crank–Nicolson scheme, and consider the following linear system:

If *η*(*t*) = *η*_{0} + *δη*(*t*), with *δη*(0) = 0, and if a Crank–Nicolson scheme is applied with a time step *δt*, the discretized equation reads

which gives the time evolution of *η*, since *δ̂η*(*δt*) ≈ *δη*(*δt*) for small *δt*. For any *t* > 0, *δ̂η*(*t*) is given by

The Borel transform of the system (A8) reads, because ℬ(*k*) = *k* for a constant function *k*,

###### Time evolution of the model

###### TLS analysis

As is well known, poles of Laplace transform of TLS solutions are eigenmodes of the system. The same holds for Borel transform: determining the poles of the Borel transform yields the complete dynamics of the system.

ZOOM is able to compute numerically the Borel transform of the TLS solution {ℬ[** δ̃η**](

*τ*) and ℬ[

**](**

*δ̃ϕ**τ*)} on the real axis

*τ*> 0. The problem of describing the dynamics of a system is thus reduced to that of determining the poles of the Borel transform of the TLS solution from its numerical values on the positive real axis.

In particular, in Eq. (A7), the poles of ℬ[** δ̃ϕ**](

*τ*) are (i) the poles of ℬ[

*δ̃ϕ*_{ins}], that is, the poles of the model without taking into account the interactions between subsystems; (ii) the poles of (1 +

**C**)

^{−1}, that is, the poles corresponding to subsystem interaction. The inverse Borel transform of Eq. (A7), obtained by an identification of simple elements, provides full dynamics of the model. The methodology consists here of fitting the Borel transform with a linear combination of sigmoid and bump functions, which are the only possible Borel transforms of linear differential equation solutions. From the characteristic times of the corresponding poles and their residue, the original function can easily be reconstructed without inverse Borel transform.

### APPENDIX B

#### Implementation of the Model under the TEF

Under the TEF, the system summarized in Fig. 1 is composed of five subsystems (five cells) listed in Table B1, together with their state variables and equations. The connections between these subsystems are represented by seven “transfer” models (13 transfer variables), listed in Table B2. Parameters are set as given in Tables B3 and B4. The role of the temperature *T*_{WV} is detailed in section 3. A nomenclature is reproduced in Table B5.

## Footnotes

*Corresponding author address:* S. Hallegatte, Center for Environment Sciences and Policy, Stanford Institute for International Studies, Encina Hall E504A, 616 Serra Street, Stanford, CA 94305. Email: hallegatte@centre-cired.fr

* Centre International de Recherche sur l'Environnement et le Développement is a joint laboratory of the Centre National de la Recherche Scientifique, l'Ecole des Hautes Etudes en Sciences Sociales, l'Ecole Nationale des Ponts-et-Chaussées, and l'Ecole Nationale du Génie Rural, des Eaux et des Forêts.

^{1}

As we will see, such a suppression requires the ability to unambiguously distinguish the feedback loop from the rest of the model.

^{2}

We assume there is no water vapor in the stratosphere.

^{3}

Note that we derive the TLS around the current state of the model, not around a reference trajectory.

^{4}

The Borel transform is similar to Laplace transform, the Borel transform variable being the inverse of the Laplace variable and homogeneous to time.

^{5}

The same type of problem occurs in GCM feedback analyses and explains why it is very difficult to carry out rigorous feedback analyzes in GCM (see Schneider et al. 1999). For example, during a CO_{2} increase simulation, cutting out WV feedback by replacing humidity by some control humidity fields entails the lost of short-term local consistency between cloud cover and humidity.

^{6}

The unit step function is only needed for the mathematical definition of the feedback function, the response to realistic perturbations can easily be derived, as we will show in the following section.

^{7}

As a comparison, the time needed by the ocean mixed layer to adjust to a flux of 1 W m^{−2} with a final temperature change of 1 K is 6.3 yr.

^{8}

This ramp can be viewed as a very stylized modeling of the global warming since 1950.

^{9}

The static influence of each pole constitutes its influence when all transients have been damped out.