3.2 Deterministic and Stochastic Processes

Learning Objectives

After completing this lecture, students should be able to:

explain the difference between deterministic and stochastic processes
interpret the role of damping in time-evolution equations
describe how stochastic forcing represents unresolved variability
formulate and implement simple deterministic and stochastic models
explain the concept of a Markov process and its implications for predictability

Concept

Climate variability arises from a combination of deterministic dynamics and stochastic processes.

We describe the time evolution of a system using a tendency equation, written in continuous form as a partial differential equation:

\[ \frac{\partial T}{\partial t} = F \]

Even though $T$ may depend on space and time in general, in this lecture we consider a simplified case in which the system depends only on time. This allows us to focus on the role of different types of forcing in time evolution.

A deterministic process evolves in a predictable way once the initial condition is given. For example, a damped system gradually relaxes toward an equilibrium state, and its future evolution is fully determined by its governing equation.

However, the climate system also contains many unresolved processes, such as small-scale turbulence and weather variability. These processes introduce randomness that cannot be explicitly predicted. To represent their effects, we use stochastic forcing.

A stochastic process includes a random component, meaning that even with identical initial conditions, the system can evolve differently. This randomness is typically modeled as noise added to the governing equation.

One important class of stochastic models is the Markov process, in which the future state depends only on the present state and not on the past history. This assumption simplifies modeling and is widely used in climate dynamics and predictability studies.

In this lecture, we first introduce a damped deterministic system, then extend it by adding stochastic forcing, and finally interpret the resulting system as a Markov (AR1) process.

Model / Equation

As discussed in the previous lecture, we consider a tendency equation in which both the state variable and forcing depend only on time $t$:

\[ \frac{\partial T}{\partial t} = F \]

In this lecture, we examine three types of forcing:

deterministic (damped) forcing
stochastic forcing
their combination, known as an auto-regressive process of order 1 (AR(1)) or a Markov process

For numerical integration, we use the forward time scheme (Euler method). The discrete form is:

\[ T^{n+1} = T^n + \Delta t \cdot F^n \]

Damped System (Deterministic Process)

We begin with a simple damped system:

\[ F = -\lambda \, T \]

where:

$T$ : state variable (e.g., temperature anomaly)
$\lambda$ : damping coefficient

The governing equation becomes:

\[ \frac{\partial T}{\partial t} = -\lambda \, T \]

The analytical solution is:

\[ T(t) = T_0 \, e^{-\lambda t} \]

This represents exponential decay toward equilibrium.

The numerical solution using the forward scheme is:

\[ T^{n+1} = T^n \Delta t \cdot \left( - \lambda \, T^n \right) \]

This system is fully determined by the governing equation and initial condition, and is therefore referred to as a deterministic process.

Stochastic Forcing (Random Process)

To represent unresolved variability, we introduce a stochastic forcing:

\[ F = \epsilon \]

where:

$\epsilon$ : random forcing (noise term)

The governing equation becomes:

\[ \frac{\partial T}{\partial t} = \epsilon \]

In discrete form:

\[ T^{n+1} = T^n + \Delta t \cdot \epsilon^n \]

where $\epsilon^n$ is a random variable (e.g., Gaussian noise with zero mean).

A system that includes random forcing is called a stochastic process, meaning that its evolution is not uniquely determined by the initial condition.

AR(1) Process (Damped + Stochastic)

The combination of damping and stochastic forcing is known as an auto-regressive process of order 1 (AR(1)), and is often interpreted as a Markov process:

\[ F = -\lambda \, T + \epsilon \]

The governing equation becomes:

\[ \frac{\partial T}{\partial t} = -\lambda \, T + \epsilon \]

The numerical form is:

\[ T^{n+1} = T^n + \Delta t \cdot \left(-\lambda T^n + \epsilon^n \right) \]

This system includes both deterministic decay and random forcing, and therefore represents a stochastic process with memory, in which the future state depends only on the current state.

Because the system includes random forcing, the exact trajectory cannot be determined analytically. Instead, we solve the system numerically and analyze its statistical properties.

Sample Code

Python scripts for this lecture are provided separately.

Initial template:

3_2-force_template.py

This script includes:

parameter settings
a forcing function force_term(T)
Euler time integration
a two-panel plot for forcing and temperature
an animation of the time evolution
automatic output of both .mp4 and .png files

Example command:

python 3_2-force_template.py

The script produces:

an animation of the time evolution (.mp4)
a final figure showing forcing and temperature (.png)

Exercise

Follow the steps below to implement and explore the damped deterministic system.

Step 1 Dampled process

Version 1: Implement the Damped Forcing

Start from the initial template and create a copy:

$ cp 3_2-force_template.py 3_2-damp_v1.py

Modify the forcing function to include damping:

def force_term(T):
    F = -0.1 * T
    return F

Run the script and examine:

the evolution of temperature
how the system approaches equilibrium

Focus on the following questions:

How does the temperature change over time?
Does the system approach a steady state?
What is the role of the damping term?

Version 2: Move the Parameter Outside the Function

Create a new version of the script:

$ cp 3_2-damp_v1.py 3_2-damp_v2.py

Define the damping parameter at the top:

lam = 0.1

Then modify the forcing function and its calling:

def force_term(T, lam):
    F = -1.0 * lam * T
    return F

F_save[0] = force_term(T, lam)

for it in range(ntime - 1):
    F = force_term(T, lam)

Run the script again and confirm that the result is unchanged.

Focus on the following questions:

Why is it better to define parameters outside the function?
How does this improve readability and flexibility of the code?

Version 3: Compare Different Damping Strengths

Create another version:

$ cp 3_2-damp_v2.py 3_2-damp_v3.py

Test at least two different values of the damping coefficient:

lam1 = 0.1
lam2 = 0.2

Then revise T, F, T_save, and F_save so that both damping cases are integrated in parallel:

# ==============================
# Time integration
# ==============================
t_save = np.zeros(ntime)
T1_save = np.zeros(ntime)
F1_save = np.zeros(ntime)
T2_save = np.zeros(ntime)
F2_save = np.zeros(ntime)

# Initial values
t = t0
T1 = T0
T2 = T0

t_save[0] = t
T1_save[0] = T1
T2_save[0] = T2
F1_save[0] = force_term(T1, lam1)
F2_save[0] = force_term(T2, lam2)

# Do-loop
for it in range(ntime - 1):
    F1 = force_term(T1, lam1)
    F2 = force_term(T2, lam2)
    t = t + dt
    T1 = T1 + dt * F1
    T2 = T2 + dt * F2

    t_save[it + 1] = t
    T1_save[it + 1] = T1
    F1_save[it + 1] = F1
    T2_save[it + 1] = T2
    F2_save[it + 1] = F2

We also need to revise the figure section so that both cases are plotted together. Replace the panel-plot part with the following:

# ==============================
# Two-panel figure
# ==============================
fig, axes = plt.subplots(2, 1, figsize=(10, 7), sharex=True)

# Group data for simpler plotting
F_save_list = [F1_save, F2_save]
T_save_list = [T1_save, T2_save]
labels = [f'lam={lam1}', f'lam={lam2}']
colors = ['C1', 'C2']
markers = ['o-', 's-']
point_markers = ['o', 's']

# Create line and point objects
linesF = []
pointsF = []
linesT = []
pointsT = []

for i in range(len(labels)):
    lineF, = axes[0].plot([], [], markers[i], color=colors[i], label=f'Forcing ({labels[i]})')
    pointF, = axes[0].plot([], [], point_markers[i], color=colors[i])
    linesF.append(lineF)
    pointsF.append(pointF)

    lineT, = axes[1].plot([], [], markers[i], color=colors[i], label=f'Temperature ({labels[i]})')
    pointT, = axes[1].plot([], [], point_markers[i], color=colors[i])
    linesT.append(lineT)
    pointsT.append(pointT)

axes[0].set_ylabel('F')
axes[0].set_title('Forcing')
axes[0].grid(True)
axes[0].legend()

axes[1].set_xlabel('Time')
axes[1].set_ylabel('T')
axes[1].set_title('Temperature')
axes[1].grid(True)
axes[1].legend()

# Axis ranges
fmin = min(np.min(F) for F in F_save_list)
fmax = max(np.max(F) for F in F_save_list)
tmin = min(np.min(T) for T in T_save_list)
tmax = max(np.max(T) for T in T_save_list)

if fmax == fmin:
    fpad = 0.5
else:
    fpad = 0.1 * (fmax - fmin)

if tmax == tmin:
    tpad = 0.5
else:
    tpad = 0.1 * (tmax - tmin)

axes[0].set_xlim(t_save[0], t_save[-1])
axes[0].set_ylim(fmin - fpad, fmax + fpad)

axes[1].set_xlim(t_save[0], t_save[-1])
axes[1].set_ylim(tmin - tpad, tmax + tpad)


# ==============================
# Animation
# ==============================
def init():
    for line in linesF + linesT:
        line.set_data([], [])
    for point in pointsF + pointsT:
        point.set_data([], [])
    return linesF + pointsF + linesT + pointsT


def update(frame):
    for i in range(len(labels)):
        linesF[i].set_data(t_save[:frame + 1], F_save_list[i][:frame + 1])
        pointsF[i].set_data([t_save[frame]], [F_save_list[i][frame]])

        linesT[i].set_data(t_save[:frame + 1], T_save_list[i][:frame + 1])
        pointsT[i].set_data([t_save[frame]], [T_save_list[i][frame]])

    return linesF + pointsF + linesT + pointsT


ani = animation.FuncAnimation(
    fig,
    update,
    frames=ntime,
    init_func=init,
    interval=300,
    blit=True
)

plt.tight_layout()

script_name = os.path.splitext(os.path.basename(__file__))[0]
outfile_mp4 = script_name + '.mp4'

writer = animation.FFMpegWriter(fps=4)
ani.save(outfile_mp4, writer=writer)
print('Animation saved as:', outfile_mp4)

outfile_png = script_name + '.png'
plt.savefig(outfile_png, dpi=150)
print('Figure saved as:', outfile_png)

plt.show()

Run the script for each case and compare the results.

Focus on the following questions:

How does the damping coefficient affect the rate of decay?
Which case reaches equilibrium faster?
How does the damping timescale relate to $\lambda$?

This experiment demonstrates how the damping parameter controls the memory and persistence of the system.

Step 2 Random Process

Version 1 Implement Random Forcing

Start from the initial template and create a copy:

$ cp 3_2-force_template.py 3_2-random_v1.py

NumPy provides several functions to generate random numbers.
One of the simplest is:

np.random.random_sample() → returns a random number between 0 and 1 (uniform distribution)

Modify the forcing function:

def force_term(T):
    F = np.random.random_sample()
    return F

Run the script and examine the results:

The forcing is always positive (between 0 and 1)
The temperature tends to increase over time due to positive forcing

To center the forcing around zero, shift the range from [0, 1] to [-0.5, 0.5]:

def force_term(T):
    F = np.random.random_sample() - 0.5
    return F

Run the script multiple times and observe how the results differ each time.

Version 2 Random Seed

Because the forcing is random, the result changes every time the script is executed.
To reproduce the same result, we can fix the random seed.

$ cp 3_2-random_v1.py 3_2-random_v2.py

Add the following line at the top of the script:

np.random.seed(100)
def force_term(T):
    F = np.random.random_sample() - 0.5
    return F

Run the script multiple times and confirm:

The result is identical each time
The randomness is now reproducible

Change the seed value to generate a different realization:

$ cp 3_2-random_v2.py 3_2-random_v2_2.py

np.random.seed(50)

Run again and compare the results.

Version 3 Common Random Distributions

There are two commonly used random distributions:

Uniform distribution → values are equally likely within a specified range
Normal (Gaussian) distribution → values cluster around a mean with a given spread

Create a new version:

$ cp 3_2-random_v2.py 3_2-random_v3.py

Add histograms of both the forcing and the response:

# ==============================
# Histogram figure
# ==============================
fig_hist, axes_hist = plt.subplots(2, 1, figsize=(6, 7))

# Histogram of forcing
axes_hist[0].hist(F_save, bins=15)
axes_hist[0].set_xlabel('F')
axes_hist[0].set_ylabel('Frequency')
axes_hist[0].set_title('Histogram of Forcing')
axes_hist[0].grid(True)

# Histogram of temperature
axes_hist[1].hist(T_save, bins=15)
axes_hist[1].set_xlabel('T')
axes_hist[1].set_ylabel('Frequency')
axes_hist[1].set_title('Histogram of Temperature')
axes_hist[1].grid(True)

plt.tight_layout()

# Save histogram figure
outfile_hist = script_name + '_hist.png'
plt.savefig(outfile_hist, dpi=150)
print('Histogram saved as:', outfile_hist)

Then run the code. The sample size is still relatively small, so increase the integration length by modifying ntime:

ntime = 100

Check the histograms of both the forcing and the response.

Uniform distribution:

np.random.seed(50)

def force_term(T):
    #F = np.random.random_sample() - 0.5
    F = np.random.uniform(-1.0, 1.0)
    return F

np.random.uniform(a, b) generates values between a and b
All values within this range are equally likely

Run the script and examine the forcing and temperature distributions.

Normal (Gaussian) distribution:

np.random.seed(50)

def force_term(T):
    #F = np.random.random_sample() - 0.5
    #F = np.random.uniform(-1.0, 1.0)
    F = np.random.normal(0.0, 1.0)
    return F

np.random.normal(mean, std) generates values with:
mean = 0.0
standard deviation = 1.0
Values near the mean are more likely than extreme values

Run the script and compare:

How does the forcing distribution differ from the uniform case?
How does the temperature distribution differ from the forcing distribution?
How does the choice of forcing distribution affect the evolution of temperature?

Step 3 AR(1) Process

The AR(1) process is a linear combination of a damped process and random forcing.

To examine the response sensitivity to the damping parameter, we will reuse the two-line plotting structure from the damped-process exercise.

$ cp 3_2-damp_v3.py 3_2-ar1_v1.py

Modify the forcing function as follows:

def force_term(T, lam):
    F = -lam * T + np.random.normal(0.0, 1.0)
    return F

Run the code and examine the results.

Because this approach generates independent random forcing for each case, the effect of the damping parameter may be difficult to isolate. The differences between the two simulations may reflect both the damping coefficient and the random forcing.

Version 2 Use the Same Random Forcing

To compare the effect of the damping coefficient more clearly, modify the code so that both simulations use the same random forcing at each time step.

$ cp 3_2-ar1_v1.py 3_2-ar1_v2.py

Set the initial condition to zero:

T0 = 0.0

Then separate the damped forcing from the random forcing. First, return only the damped term in the function:

def force_term(T, lam):
    F = -lam * T
    return F

Next, add the same random forcing Fr to both cases inside the time-integration loop:

for it in range(ntime - 1):
    Fr = np.random.normal(0.0, 1.0)
    F1 = force_term(T1, lam1) + Fr
    F2 = force_term(T2, lam2) + Fr

Run the code again and compare the results.

Now the two simulations are driven by the same random forcing, so the differences between them more clearly reflect the impact of the damping coefficient.

Version 3 Explore Parameter Sensitivity

Try several combinations of the damping coefficient and random forcing strength.

For example, modify:

lam1 and lam2
the standard deviation of the random forcing in np.random.normal(0.0, 1.0)

Focus on the following questions:

How does the damping coefficient affect the persistence of the response?
How does stronger random forcing affect the variability of the system?
Which case shows longer memory?
How does the balance between damping and random forcing control the overall behavior?

This experiment demonstrates that the AR(1) process represents a system with both memory and random variability.

Key Points

The tendency equation provides a general framework for describing time evolution: $$ \frac{\partial T}{\partial t} = F $$
Deterministic processes are fully determined by the governing equation and initial condition, while stochastic processes include random forcing and are inherently unpredictable in detail.
A damped system ($F = -\lambda T$) exhibits exponential decay, where the damping coefficient $\lambda$ controls the rate of adjustment and memory of the system.
Stochastic forcing represents unresolved variability and introduces randomness into the system, leading to different realizations even with identical initial conditions.
The choice of random distribution (uniform vs normal) affects the statistical properties of the forcing and the system response.
The AR(1) process combines damping and stochastic forcing: $$ F = -\lambda T + \epsilon $$ and represents a system with both memory and variability.
In the AR(1) system, the balance between damping and stochastic forcing determines:
persistence (memory)
variability (spread of the response)
Using identical random forcing across simulations is essential for isolating parameter sensitivity.

Homework

Implement the AR(1) model using two different damping coefficients (e.g., $\lambda_1$ and $\lambda_2$), and use the same random forcing for both cases.

Plot the time evolution of temperature for both cases.
Generate histograms of temperature.
Briefly answer (3–5 sentences):

Which case shows stronger persistence (memory)?
How does the damping coefficient affect variability?

Submit:

one figure (time series and histogram)
a short written explanation