Gradient aware harmonisation (Getting Started)¶

Here we introduce a method for harmonising two timeseries. This part may be more unusual or unfamiliar to people used to working with arrays, so it serves as an introduction into some of the concepts used in this package.

Imports¶

In [1]:

# Imports
import matplotlib.pyplot as plt
import numpy as np

from gradient_aware_harmonisation import harmonise
from gradient_aware_harmonisation.helpers import plotting
from gradient_aware_harmonisation.utils import (
    Timeseries,
    cosine_decay,
    polynomial_decay,
)

# Imports import matplotlib.pyplot as plt import numpy as np from gradient_aware_harmonisation import harmonise from gradient_aware_harmonisation.helpers import plotting from gradient_aware_harmonisation.utils import ( Timeseries, cosine_decay, polynomial_decay, )

Toy Example 1: Artificial data¶

Create some data¶

In [2]:

# create some data
harmonisation_time = 3.0

x1 = np.arange(-2, 3.0, 0.1)
y1 = -16 * x1

x2 = np.arange(2, 10, 0.1)
y2 = 0.5 * x2 + x2**3

target_timeseries = Timeseries(time_axis=x1, values=y1)
harmonisee_timeseries = Timeseries(time_axis=x2, values=y2)

# create some data harmonisation_time = 3.0 x1 = np.arange(-2, 3.0, 0.1) y1 = -16 * x1 x2 = np.arange(2, 10, 0.1) y2 = 0.5 * x2 + x2**3 target_timeseries = Timeseries(time_axis=x1, values=y1) harmonisee_timeseries = Timeseries(time_axis=x2, values=y2)

Harmonise functions using harmonise¶

Inspect harmoniser function¶

In [3]:

help(harmonise)

harmonised_timeseries = harmonise(
    target_timeseries=target_timeseries,
    harmonisee_timeseries=harmonisee_timeseries,
    harmonisation_time=harmonisation_time,
    convergence_time=None,
    decay_method="cosine",
)

help(harmonise) harmonised_timeseries = harmonise( target_timeseries=target_timeseries, harmonisee_timeseries=harmonisee_timeseries, harmonisation_time=harmonisation_time, convergence_time=None, decay_method="cosine", )

Help on function harmonise in module gradient_aware_harmonisation:

harmonise(target_timeseries: gradient_aware_harmonisation.utils.Timeseries, harmonisee_timeseries: gradient_aware_harmonisation.utils.Timeseries, harmonisation_time: Union[int, float], convergence_time: Union[int, float, NoneType], interpolation_target: str = 'original', decay_method: str = 'cosine', **kwargs: Any) -> gradient_aware_harmonisation.utils.Timeseries
    Harmonise two timeseries
    
    When we say harmonise, we mean make it
    such that the harmonisee matches with the target at some
    specified time point (called harmonisation time)
    
    Parameters
    ----------
    target_timeseries
        Target timeseries (i.e. what we harmonise to)
    
    harmonisee_timeseries
        Harmonisee timeseries (i.e. the timeseries we want to harmonise)
    
    harmonisation_time
        Time point at which harmonisee should be matched to the target
    
    convergence_time
        Time point at which the harmonised data
        should converge towards the prediced data.
    
    interpolation_target
        Target to which the harmonised timeseries should converge.
    
        If original, we converge back to harmonisee.
        If bias-corrected, we converge back to harmonissee
        having applied a basic constant offset bias correction
        (see the docs for further info TODO put a cross link to a notebook).
    
    decay_method
        Decay function used to decay weights
        when interpolating between the target and our harmonisation target.
        If 'polynomial' is used an additional argument 'pow' to specify the power
        is required (should be => 1.)
    
    **kwargs
        keyword arguments passed to `make_interp_spline` or 'polynomial_decay'
    
    Returns
    -------
    harmonised_timeseries :
        timeseries of harmonised data set
    
    Raises
    ------
    ValueError
        interpolation_target must be either 'original' or 'bias_corrected'

Plot results¶

In [4]:

plotting(
    harmonisee_timeseries,
    target_timeseries,
    harmonised_timeseries,
    harmonisation_time,
    convergence_time=None,
)

plotting( harmonisee_timeseries, target_timeseries, harmonised_timeseries, harmonisation_time, convergence_time=None, )

Introduce a convergence_time at $x=8.$¶

In [5]:

convergence_time = 8.0

harmonised_timeseries = harmonise(
    target_timeseries=target_timeseries,
    harmonisee_timeseries=harmonisee_timeseries,
    harmonisation_time=harmonisation_time,
    convergence_time=convergence_time,
    decay_method="cosine",
)

plotting(
    harmonisee_timeseries,
    target_timeseries,
    harmonised_timeseries,
    harmonisation_time,
    convergence_time,
)

convergence_time = 8.0 harmonised_timeseries = harmonise( target_timeseries=target_timeseries, harmonisee_timeseries=harmonisee_timeseries, harmonisation_time=harmonisation_time, convergence_time=convergence_time, decay_method="cosine", ) plotting( harmonisee_timeseries, target_timeseries, harmonised_timeseries, harmonisation_time, convergence_time, )

Change weight decay_method¶

Inspect different decay variants¶

In [6]:

cos_decay = cosine_decay(decay_steps=50)
poly_decay = lambda pow: polynomial_decay(decay_steps=50, pow=pow)  # noqa: E731

plt.figure(figsize=(4, 3))
plt.plot(cos_decay, label="cosine")
for p in [1.0, 2.0, 3.0]:
    plt.plot(poly_decay(p), label=f"polynomial, pow={p}")
plt.legend(handlelength=0.3, fontsize="small", frameon=False)
plt.title("Different weight decay methods")
plt.show()

cos_decay = cosine_decay(decay_steps=50) poly_decay = lambda pow: polynomial_decay(decay_steps=50, pow=pow) # noqa: E731 plt.figure(figsize=(4, 3)) plt.plot(cos_decay, label="cosine") for p in [1.0, 2.0, 3.0]: plt.plot(poly_decay(p), label=f"polynomial, pow={p}") plt.legend(handlelength=0.3, fontsize="small", frameon=False) plt.title("Different weight decay methods") plt.show()

Inspect convergence results using different decay methods¶

In [7]:

def sensitivity_weight_decay(decay_method, **kwargs):  # noqa: D103
    harmonised_timeseries = harmonise(
        target_timeseries=target_timeseries,
        harmonisee_timeseries=harmonisee_timeseries,
        harmonisation_time=harmonisation_time,
        convergence_time=None,
        decay_method=decay_method,
        **kwargs,
    )

    plotting(
        harmonisee_timeseries,
        target_timeseries,
        harmonised_timeseries,
        harmonisation_time,
        convergence_time=None,
    )

def sensitivity_weight_decay(decay_method, **kwargs): # noqa: D103 harmonised_timeseries = harmonise( target_timeseries=target_timeseries, harmonisee_timeseries=harmonisee_timeseries, harmonisation_time=harmonisation_time, convergence_time=None, decay_method=decay_method, **kwargs, ) plotting( harmonisee_timeseries, target_timeseries, harmonised_timeseries, harmonisation_time, convergence_time=None, )

Cosine decay¶

In [8]:

sensitivity_weight_decay(decay_method="cosine")

sensitivity_weight_decay(decay_method="cosine")

Polynomial decay, power = 1. (i.e., linear decay)¶

In [9]:

sensitivity_weight_decay(decay_method="polynomial", pow=1.0)

sensitivity_weight_decay(decay_method="polynomial", pow=1.0)

Polynomial decay, power = 2.¶

In [10]:

sensitivity_weight_decay(decay_method="polynomial", pow=2.0)

sensitivity_weight_decay(decay_method="polynomial", pow=2.0)

Polynomial decay, power = 3.¶

In [11]:

sensitivity_weight_decay(decay_method="polynomial", pow=3.0)

sensitivity_weight_decay(decay_method="polynomial", pow=3.0)

Toy Example 2: Use timeseries data¶

Read data¶

In [12]:

# select harmonization time point
harmonisation_time = 2004
convergence_time = 2006

# get timeseries
target_timeseries = Timeseries(
    time_axis=np.array([2001, 2002, 2003, 2004]),
    values=np.array([371.77, 373.72, 376.33, 378.83]),
)
harmonisee_timeseries = Timeseries(
    time_axis=np.array([2003, 2004, 2005, 2006, 2007]),
    values=np.array([375.56, 376.28, 378.83, 381.20, 382.55]),
)

# select harmonization time point harmonisation_time = 2004 convergence_time = 2006 # get timeseries target_timeseries = Timeseries( time_axis=np.array([2001, 2002, 2003, 2004]), values=np.array([371.77, 373.72, 376.33, 378.83]), ) harmonisee_timeseries = Timeseries( time_axis=np.array([2003, 2004, 2005, 2006, 2007]), values=np.array([375.56, 376.28, 378.83, 381.20, 382.55]), )

Run harmonise for different settings¶

In [13]:

# harmonise timeseries at t0
harmonised_timeseries = harmonise(
    target_timeseries,
    harmonisee_timeseries,
    harmonisation_time=harmonisation_time,
    convergence_time=None,
)

# harmonise timeseries at t0 and assure convergence at t1 (converge_t)
harmonised_timeseries2 = harmonise(
    target_timeseries,
    harmonisee_timeseries,
    harmonisation_time=harmonisation_time,
    convergence_time=convergence_time,
)

# harmonise timeseries at t0 and assure convergence at t1 (converge_t)
harmonised_timeseries3 = harmonise(
    target_timeseries,
    harmonisee_timeseries,
    harmonisation_time=harmonisation_time,
    convergence_time=convergence_time,
    interpolation_target="bias_corrected",
)

# harmonise timeseries at t0 harmonised_timeseries = harmonise( target_timeseries, harmonisee_timeseries, harmonisation_time=harmonisation_time, convergence_time=None, ) # harmonise timeseries at t0 and assure convergence at t1 (converge_t) harmonised_timeseries2 = harmonise( target_timeseries, harmonisee_timeseries, harmonisation_time=harmonisation_time, convergence_time=convergence_time, ) # harmonise timeseries at t0 and assure convergence at t1 (converge_t) harmonised_timeseries3 = harmonise( target_timeseries, harmonisee_timeseries, harmonisation_time=harmonisation_time, convergence_time=convergence_time, interpolation_target="bias_corrected", )

Plot results¶

In [14]:

plotting(
    harmonisee_timeseries=harmonisee_timeseries,
    target_timeseries=target_timeseries,
    interpolated_timeseries=harmonised_timeseries,
    harmonisation_time=harmonisation_time,
    convergence_time=None,
)

plotting(
    harmonisee_timeseries,
    target_timeseries,
    harmonised_timeseries2,
    harmonisation_time,
    convergence_time,
)

plotting(
    harmonisee_timeseries,
    target_timeseries,
    harmonised_timeseries3,
    harmonisation_time,
    convergence_time,
)

plotting( harmonisee_timeseries=harmonisee_timeseries, target_timeseries=target_timeseries, interpolated_timeseries=harmonised_timeseries, harmonisation_time=harmonisation_time, convergence_time=None, ) plotting( harmonisee_timeseries, target_timeseries, harmonised_timeseries2, harmonisation_time, convergence_time, ) plotting( harmonisee_timeseries, target_timeseries, harmonised_timeseries3, harmonisation_time, convergence_time, )