Timescales of Anthropogenic Perturbation to Earth System Processes

4.1 General characterization of age and transit time distributions

Consider an open system, potentially containing multiple interconnected compartments, that receives and releases material through its boundaries. This material is characterized by its mass, and we will consider the total mass stored in the system, in the input flux, and in the output flux. We will also consider these masses as composed of small units that could be interpreted as individual atoms or molecules, which, for generality, I call particles from now on. In addition, I consider two different representations of time, a time dimension T that serves as a reference frame for the observer and can be expressed as a calendar time; and as time intervals between events occurring on T (Figure 1, Table 1). When particles enter the system, the event of their time of entry or injection is recorded as t_i, and the time when they exit as t_e. Thus, the time interval $\tau =t_e-t_i$ is defined as the transit time of individual particles traveling through the system. Observing the system is also an event, and it occurs at t. Thus, we are interested in the time interval between the time of entry and the time of observation of a particle inside the system, $a=t-t_i$, and we define this interval as the age of the particle. The interval between observing the particle in the system and its time of exit, $s=t_e-t$, is defined as the survival time, which has also been called residual or remaining life expectancy (Nauman, 1969, 2008). Although the survival time has not been addressed in detail in most previous publications (however see Calabrese and Porporato, 2015), it emerges from the same conceptual basis as age and transit time. These three metrics are related by the expression $\tau =a+s$ (Figure 1).

The theory of timescale characteristics in Earth system processes is based on the macroscopic characterization of the time intervals just defined. These intervals are treated as random variables for which probability distributions and their moments can be computed. The traditional theory for systems in equilibrium is only based on the characterization of the properties of these intervals, and it assumes that the time when particles enter is irrelevant. For systems out of equilibrium, the theory considers both the reference time frame where events occur and the intervals between events (Rasmussen et al., 2016).

In the general theory for systems out of equilibrium, we are interested in describing the timescale characteristics of the mass stored in the system at any given time M(t). Because we assume that the total mass is composed of particles that have entered the system in the past and have a certain age at the time of observation, we can express the total mass as a cumulative distribution function M(a,t) over the variable a for fixed values of t (i.e., a marginal distribution). We assume that as we accumulate mass over increasing values of age, the function has a well-defined limit

Equivalently, the total mass in the system at t can be expressed in terms of a pdf of ages as $\psi (a,t)$

where the prime symbol in a’ indicates the dummy variable of integration, and ${\int }_0^{\mathrm{\infty }}\psi (a,t)\text{d}a=1$ for all values of t. Differentiating this equation with respect to a and rearranging terms yields

Thus, the pdf of ages of mass stored in the system can be interpreted as the relative sensitivity of the cumulative mass function with respect to ages.

A similar characterization can be made for the fluxes in or out of the system with respect to the variable transit time $\tau$. However, for systems out of equilibrium, it is necessary to differentiate between the transit time of particles arriving at t that will have a future transit time, versus the particles being released out of the system at time t that had a past transit time. We refer to these different forms of transit time as the backward and the forward transit time, respectively. For simplicity of exposition, we will focus for the moment on the backward transit time and assume that the flux F(t) represents the output flux from the system at time t (in section 4.4, the theory to the forward case is expanded). In analogy to mass storage, this flux can be characterized by a cumulative distribution function with an upper limit

which can also be expressed in terms of a pdf

with ${\int }_0^{\mathrm{\infty }}\phi (\tau ,t)\text{d}\tau =1$. This equation can be manipulated to express the transit time pdf as

Similarly, a cumulative distribution function and a pdf for the survival time can be defined as for age and transit time. However, they are not presented here due to space limitations and because they are not relevant to the theory presented below.

4.2 Relationships between age and transit time distributions

The age $\psi (a,t)$ and the transit time $\phi (\tau ,t)$ pdfs are of particular interest in many scientific disciplines, and most efforts shall concentrate on obtaining estimates of these functions. For the special case of systems in equilibrium, these two functions are related in a closed mathematical form (Bolin and Rodhe, 1973; Björkström, 1978). For systems out of equilibrium, their relationship is less straightforward (however, see O’Neill et al. (1994) for a closed-form relationship between the means of these distributions), and for this reason their relation is better studied using the ratios

These are different forms of the so-called SAS functions used in watershed hydrology (Benettin et al., 2022), and they can be used to infer the proportion of the storage of mass at different ages that is selected for release in the output flux. Equation (8) in particular, can be interpreted as the probability that the mass of a certain age will be released out of the system at a particular time, or as the removal rate in units of inverse time. Using equations (3) and (6), and the fact that $\partial \tau =\partial a$, we can express equation (8) as

which can be interpreted as the change in the cumulative transit time distribution at time t generated by a change in the cumulative age distribution at t.

The particular shapes of the functions of equations (7) to (9) at each value of t provide very useful information about how the output flux “selects” (as in statistical selection) mass across the range of available ages to be released out of the system. Thus, three type of shapes are possible: (a) most of the mass that appears in the output flux is composed of old age, (b) the mass in the output flux is uniformly selected across the entire range of available ages of stored mass, and (c) the mass in the output flux is preferentially selected from young ages of mass stored (Figure 4). These three cases generalize the three cases for the means described by Bolin and Rodhe (1973) for systems in equilibrium (Figure 2), and offer a powerful approach to study overall characteristics of system dynamics.

Figure 4

4.3 General solution of the time-dependent age distribution

So far, nothing has been said about how to obtain the different forms of age and transit time distributions. Working independently, researchers in hydrology (Botter et al., 2011; Harman, 2015; van der Velde et al., 2012) and in carbon cycle (Metzler et al., 2018; Rasmussen et al., 2016) have developed methods to obtain age and transit time distributions for systems out of equilibrium. The approach in hydrology uses time-series observations of precipitation and run-off and makes no assumptions on the internal structure of the system. The approach parameterizes forms of the SAS functions to obtain solutions of transit time pdfs that match the observed data. In contrast, the approach in carbon cycle science relies on an available model expressed as compartmental dynamical system that makes explicit assumptions about the structure of the system as a set of compartments that exchange mass among each other, potentially including nonlinear functions that determine rates of mass exchange and release.

Despite differences in these two approaches, the new methods to obtain time-dependent age and transit time distributions have one important aspect in common, and it is the representation of the problem in the form of the McKendrick–von Foester equation (Keyfitz and Keyfitz, 1997). In mathematical biology and demography, the McKendrick–von Foester equation generalizes problems in which the size of a population (number of individuals, expressed as a function of age and time $P(a,t)$) changes with respect to a and t simultaneously as

where $\mu (a,t)$ is the rate of removal (mortality in population studies) with respect to age a at time t. This equation expresses the fact that the rate of change of the population with respect to time at age a exactly balances the rate of change with respect to age after accounting for the rate of removals.

For mass balance problems, the McKendrick–von Foester equation can be expressed as

where we interpret $\omega (a,t)$ as the rate at which mass is selected from the total mass across the range of ages. To obtain a particular solution for this equation, it is necessary to include an initial distribution as initial condition

representing the distribution of ages of the population at $t=0$; and the boundary condition

expressing the fact that all flux inputs $F_i(t)$ (newborns in population studies) have, by definition, an age of zero.

There are a number of mathematical methods available to solve equation (11) (Keyfitz and Keyfitz, 1997), which yield an analytical or numerical expression for $m(a,t)$. Nevertheless, an expression for $\omega (a,t)$ must be provided, which for a compartmental dynamical system is straightforward to specify because it emerges from the mass balance equations out of the compartments; i.e., the term $-\omega (a,t)m(a,t)$ emerges directly from the compartmental matrix (Metzler et al., 2018; Rasmussen et al., 2016). When a model is not available, as in the case of watershed hydrology problems, the proposed approach is to use an empirically parameterized SAS function to specify how the output flux selects from the available ages in storage (Benettin et al., 2022). Readers interested in the application of these methods can obtain additional details and access to software in Metzler et al. (2018), Benettin and Bertuzzo (2018), and Harman and Xu Fei (2024).

4.4 Forward and backward transit times

Once a solution for $m(a,t)$ is obtained, and given that M(t) is known, it is then possible to calculate the time-dependent age pdf $\psi (a,t)$ and the time-dependent cumulative distribution $M(a,t)$ using the equations presented above. To obtain the time-dependent forward and backward transit time distributions, we take advantage of the complete knowledge of the time evolution of the distribution of ages, and the fact that time-dependent input $F_i(t)$ and output (exit) $F_e(t)$ fluxes are also known.

We now define the backward transit time distribution function $p_{\text{BTT}}(\tau ,t)$ as a (exit) flux-weighted pdf of transit times. From equation (8), it can be obtained as

Similarly, we define the forward transit time distribution function $p_{\text{FTT}}(\tau ,t)$ as a (input) flux-weighted pdf of transit times

which implies that we have knowledge of the future dynamics of the inputs up to an arbitrary long time $t+\tau$.

Notice that for a fixed value of t in a non-equilibrium system, $p_{\text{BTT}}(\tau ,t)\ne p_{\text{FTT}}(\tau ,t)$. This is because input and output fluxes are not necessarily the same for a given t, and the functions $\omega$ and m change over t. Conversely, for systems in equilibrium where $F_i=F_e$, and $\omega$ and m are unique, it is well known that $p_{\text{BTT}}(\tau ,t)=p_{\text{FTT}}(\tau ,t)$ (Niemi, 1977). In addition, if we focus on the behavior of one single particle injected at a fixed time t_i and exiting at a fixed time t_e, we obtain the equality (Metzler et al., 2018)

In words, the transit time for an individual particle is the same independently if one looks at the time it enters or at the time it leaves the system. This equality not only applies to one single particle, but also to the set of particles that simultaneously arrive at t_i and exit together at t_e, which could be a very small subset from all particles in the input and output fluxes. However, this behavior for one individual particle or a small subset of particles does not scale up for the total fluxes because all particles entering at one single time leave at different future times, and the flux of all particles exiting at one single time entered at different past times.

In addition to these transit time distributions, it is also possible to obtain time-dependent survival time distributions $m(s,t)$. They provide additional insights about the future behavior of mass stored in a system. The survival time distributions are at the same level of importance as age and transit time distributions and complement the information provided by them (Figure 1). These three distributions are mathematically related and therefore connected via the SAS functions. Thus, one can obtain survival time distributions from knowledge of the other complementary distributions (Calabrese and Porporato, 2015).

5.1 The transit time of anthropogenic CO₂

The question of how long would it take to remove from the atmosphere CO₂ of fossil-fuel origin? has been a main source of motivation for the development of a theory of timescales of material cycling in natural systems, particularly for the early founding figures of Tellus. For instance, Rodhe (1991, p. 4) in his memoir on Bert Bolin, writes “After Bert Bolin finished his Ph.D. in 1956, C.-G. Rossby advised him to move into another research area: atmospheric chemistry and questions specifically related to the residence time of pollutants in the atmosphere. That is when his interest in the carbon cycle started.” Therefore, I think it is fair to assert that most of the early theoretical development on the topic of timescales were motivated by the question of the timescale characteristics of anthropogenic CO₂. However, this has been proven to be a difficult problem.

Carbon in fossil fuels has been isolated from other components of the Earth system for very long times. Once this carbon is burned by humans, it enters the active part of the global carbon cycle. It enters the atmosphere as CO₂ where it remains stored for some time, and eventually, this carbon is transferred to the terrestrial biosphere, or to the oceans where it would remain stored for some period of time. This anthropogenic carbon can return to the atmosphere, either by the respiratory activity of organisms on land and oceans, or by inorganic ocean-atmosphere exchange. Therefore, the question of how long anthropogenic carbon remains in the atmosphere could be studied from two different points of view. From the point of view of the isolated atmospheric compartment and only considering the initial transfer, from burning of fossil fuels until CO₂ uptake by land or oceans; or from the point of view of a larger interconnected system including the atmosphere, land, and oceans, with a final receiver geological reservoir (Archer and Brovkin, 2008; Archer et al., 2009). Because one can study this problem from these two particular points of view, a lot of confusion and debates have existed regarding this question (Gaffin et al., 1995; O’Neill et al., 1994; Tans, 1997).

In addition to these difficulties about how to define the boundaries of the system of interest, there are additional difficulties related to the nonlinearity of carbon exchange between the atmosphere with oceans and land, as well as issues of time-dependencies of other external forcing factors. For these reasons, the problem of the time to remove carbon of fossil-fuel origin must be addressed under a broader conceptual framework, one that can deal with a general characterization of carbon compartments, and that can address explicitly nonlinearities and time-dependencies.

A reasonable approach to address this problem is to define an open compartmental system that has as its main input flux carbon of geological origin entering the atmosphere, and an output flux of carbon returning to a geological reservoir such as ocean sediments. Carbon is transferred among the main active compartments of atmosphere, land, and oceans, with equations defining the main transfer fluxes including potential feedback among reservoirs that lead to systems of nonlinear differential equations (Rodhe and Björkström, 1979). Rates of carbon transfer among reservoirs could also be affected by changes in climate, such as the rates of photosynthesis and respiration on land and oceans, or the rates of CO₂ exchange between the atmosphere and the ocean (Kaufhold et al., 2025). Many of these processes are already represented in Earth system models (ESMs) of different levels of complexity, from reduced complexity models and simple emulators to high-resolution ESMs. Many of these models are currently under active development, improving their representation of multiple processes with their nonlinearities and feedback. More importantly, most of these models can be translated to a mathematical representation according to the formalism of nonlinear non-autonomous compartmental system (Metzler et al., 2020; Sierra et al., 2018). Then, the question of how long would it take to remove carbon from fossil fuel origin from the active carbon cycle can be addressed as a problem of finding the forward or backward transit time of a nonlinear non-autonomous compartmental system (Metzler et al., 2018), linking inputs from the burning of fossil fuels to outputs to geological storage.

Under this theoretical framework, it would be possible to learn about differences in forward transit times between anthropogenic carbon emitted in the past, versus emissions in the future under different scenarios. It would be possible to learn about ages and transit times of individual compartments, or to obtain the survival time of the current burden of anthropogenic CO₂. It would also be possible to compute atmospheric warming effects attributed to emissions occurring at different times, expanding the current conceptual framework of static global warming potentials (Joos et al., 2013; Lashof and Ahuja, 1990; Rodhe, 1990). It would also be possible to learn how human activities perturb the global carbon cycle across a wide range of timescales and how these perturbations are exacerbated over time as emissions continue, or mitigated as emission policies are implemented. This type of analysis could be readily implemented using simple models such as FaIR (Leach et al., 2021) or BernSCM (Strassmann and Joos, 2018) combined with existing tools for timescales analysis, such as those developed in hydrology (Harman and Xu Fei, 2024) or for compartmental systems of the carbon cycle (Metzler et al., 2018). Overall, the implementation of such a computational framework would enable a set of new scientific questions and policy analyses.

5.2 Timescales of anthropogenic perturbations to Earth system processes

The previous example demonstrates that a general mathematical framework is essential to address a complex problem that requires a careful definition of system boundaries, a specification of nonlinearities in material exchange among reservoirs, and a clear specification of time-dependencies in forcing that alter rates of material exchange. It also shows that for such complex problems, one must also be specific about the type of questions being asked. Because the forward and backward transit times are different for systems out of equilibrium, one must be careful in either posing questions on the future dynamics of current material inputs or the past dynamics of current outputs. Given the extent and magnitude of human modification of material cycles in the Earth system, it is reasonable to expect that the age and transit time distributions for many material cycles have been drastically modified over the past one-and-a-half century. Therefore, scientific questions about the timescales of material cycles should move from attempts to identify single values of a timescale metric to new questions about the magnitude and direction of their perturbation.

Two areas with potential for new scientific developments regarding the application of timescale metrics are: the timescales of reactive nitrogen in the environment, and the modification of the timescales of methane in the atmosphere. Although there are many other areas with large potential for applications, I present here these two as clear examples where progress can be achieved within the next years.

Although the atmosphere stores large amounts of nitrogen (N), this N is not easily available for organisms on land and oceans because it is present as inert N₂ gas, with a triple covalent bond between the two N atoms that is very difficult to break. Since the beginning of the industrial revolution, humans have increasingly added reactive nitrogen to the environment by burning fossil fuels, by promoting the cultivation of species that can fix N, and by artificially producing reactive N through the Haber-Bosch process to produce fertilizers (Galloway et al., 2021). This reactive nitrogen is cycled among terrestrial ecosystems, the atmosphere, and aquatic ecosystems, in a chain of biogeochemical reactions collectively known as the nitrogen cascade (Fowler et al., 2013; Galloway et al., 2003). During the time reactive N stays in these reservoirs, it produces a number of negative environmental impacts, which include: contribution to the greenhouse effect by N₂O, increases in the production of tropospheric ozone that acts as air pollutant, acidification of precipitation, eutrophication of aquatic ecosystems, and changes in the biodiversity of terrestrial and aquatic ecosystems, among other impacts (Galloway et al., 2003). Eventually, reactive nitrogen is denitrified and returns to the atmosphere as N₂.

The time that anthropogenic reactive nitrogen spends in the natural environment is currently unknown. Although there are broad estimates of the timescales of reactive N in individual reservoirs, such as atmospheric N₂O (Prather et al., 2015) or in global soils and oceans (Fowler et al., 2013), there are no estimates linking these different reservoirs through the pathways in which reactive N moves across the Earth system. It is likely that through nonlinear interactions with the carbon and the phosphorous cycle (Gruber and Galloway, 2008), the forward transit time of reactive nitrogen will increase, spending more time in different reservoirs and exacerbating a number of environmental problems. As more reactive nitrogen is available for net primary production, more nitrogen is stored in organic compounds and could be cycled between ecosystems and food production systems, increasing its overall transit time. I believe there is an urgent need to address this problem so better policies for environmental management of N can be developed, taking advantage of the new theoretical and computational developments in the calculation of transit times, integrating new datasets and models of the global nitrogen cycle.

Another important anthropogenic perturbation is the growth rate of methane (CH₄) in the atmosphere and its “lifetime.” As opposed to CO₂, the trend of atmospheric CH₄ mole fraction has not increased steadily with emissions, but it has shown periods of rapid growth and periods of relatively constant values (Dlugokencky et al., 2011; Saunois et al., 2025; Turner et al., 2019). Although there are large uncertainties regarding the rate of increase of different emission sources, there seems to be good agreement that the main sink of CH₄, its oxidation with OH, has changed over time (Folberth et al., 2025; Lelieveld et al., 1998; Stevenson et al., 2020; Turner et al., 2019). The availability of OH has a strong influence on the so-called lifetime of methane, with a number of studies showing that as OH increases in the atmosphere, the reaction rate with CH₄ speeds up and therefore CH₄ spends less time in the atmosphere. However, the chemistry of CH₄ in the troposphere is relatively complex, with multiple interactions with other gases such as O₃, CO, NO_X, among others, which have a nonlinear influence on different reaction rates and therefore affect CH₄ lifetime (Lelieveld et al., 1998; Prather, 2007). This lifetime of CH₄ has been quantified as the ratio of the total atmospheric mass to the total sink flux, which is equivalent to the turnover time. It is well-established that this estimate of turnover time is not equivalent to the concepts of age or transit time for nonlinear non-equilibrium systems (Bolin and Rodhe, 1973; Lewis and Nir, 1978; Metzler et al., 2018; Rasmussen et al., 2016; Sierra et al., 2017), and therefore most of the literature on CH₄ lifetime has been performed with methods that do not correctly inform about the time that CH₄ molecules spend in an atmosphere out of equilibrium. This is problematic because CH₄, as one of the major greenhouse gases, is very relevant for the development of policies for climate change mitigation. In the policy arena, the concepts of global warming potentials and carbon dioxide equivalents are used extensively to design climate change mitigation actions, and these concepts rely on estimates of CH₄ lifetime for their computations. Therefore, it is very important to improve estimates of CH₄ lifetime, ideally as a CH₄ transit time or survival time distribution that changes over time, using the appropriate set of concepts and formulas for nonlinear non-equilibrium systems.

Fortunately, estimates of CH₄ transit time can be readily obtained using the recent set of methods discussed here applied to conceptual models (Prather, 1996, 2007) or simulation output (Stevenson et al., 2020). For instance, Prather (2007) used a simple nonlinear model with the set of chemical reactions linking CH₄ with CO and OH, studying the eigenvector and eigenvalue decomposition of the linearized system of equations at equilibrium. This simple model, or other more complex versions, can be used to analyze steady-state transit time distributions of CH₄ for pre-industrial conditions using the methods developed by Metzler and Sierra (2018), and performing perturbation analyses to explore how the perturbation of other chemical species such as CO or OH change the equilibrium transit times of CH₄. Furthermore, results from large-scale model inter-comparison projects (e.g. Stevenson et al., 2020) can be used to obtain time-dependent (backward) transit time distributions of CH₄. This can be done by translating the methods from hydrology that rely on time series of inputs (CH₄ emission sources) and outputs (CH₄ sinks) to obtain parameterized SAS functions and the corresponding estimates of transit times (Harman, 2015; Harman and Xu Fei, 2024). Alternatively, one can reconstruct a compartmental model from the numerical output of the models (Metzler et al., 2020) and subsequently obtain backward transit times using methods for non-autonomous compartmental systems (Metzler et al., 2018; Rasmussen et al., 2016).