PS E:\LLMS\Fine-tuning> & C:/Users/hemanthk.LAP53-FJS.000/AppData/Local/Programs/Python/Python310/python.exe e:/LLMS/Fine-tuning/2024notes/timeseries/dtw.py
Standard DTW Distance: 8.0
Path: [[5 5 5 5 5 5 5]
 [0 1 2 3 4 5 6]]
FastDTW Distance: 0.0
Path: [[0 1 1 2 3 4 5]
 [0 1 2 3 4 5 6]]
Traceback (most recent call last):
  File "e:\LLMS\Fine-tuning\2024notes\timeseries\dtw.py", line 378, in <module>
    distance_sub, path_sub = subsequence_dtw(series_a, series_b)
  File "e:\LLMS\Fine-tuning\2024notes\timeseries\dtw.py", line 203, in subsequence_dtw
    raise ValueError("Subsequence length must be less than or equal to the main series length.")
ValueError: Subsequence length must be less than or equal to the main series length.
PS E:\LLMS\Fine-tuning> 

import numpy as np
import numba
from numba import njit, prange
from typing import Tuple, List, Optional
import warnings


@njit
def backtrace(trace: np.ndarray) -> np.ndarray:
    """
    Backtraces the optimal path from the DTW trace matrix.

    Parameters:
        trace (np.ndarray): Trace matrix indicating the path.

    Returns:
        np.ndarray: Optimal warping path as a 2xN array.
    """
    i, j = trace.shape[0] - 1, trace.shape[1] - 1
    path = []

    while i > 0 and j > 0:
        path.append((i - 1, j - 1))
        tb = trace[i, j]

        if tb == 0:
            i -= 1
            j -= 1
        elif tb == 1:
            i -= 1
        elif tb == 2:
            j -= 1
        else:
            # Unexpected value in trace matrix
            raise ValueError(f"Unexpected trace value {tb} at position ({i}, {j})")

    path = path[::-1]
    return np.array(path).T  # Shape: (2, N)


@njit(parallel=True)
def dtw_cpu(x: np.ndarray) -> np.ndarray:
    """
    Computes the standard DTW alignment path on CPU.

    Parameters:
        x (np.ndarray): Cost matrix.

    Returns:
        np.ndarray: Optimal warping path.
    """
    N, M = x.shape
    cost = np.full((N + 1, M + 1), np.inf, dtype=np.float32)
    trace = -np.ones((N + 1, M + 1), dtype=np.int32)

    cost[0, 0] = 0

    for i in prange(1, N + 1):
        for j in range(1, M + 1):
            c0 = cost[i - 1, j - 1]
            c1 = cost[i - 1, j]
            c2 = cost[i, j - 1]

            if c0 < c1 and c0 < c2:
                cost[i, j] = x[i - 1, j - 1] + c0
                trace[i, j] = 0
            elif c1 < c0 and c1 < c2:
                cost[i, j] = x[i - 1, j - 1] + c1
                trace[i, j] = 1
            else:
                cost[i, j] = x[i - 1, j - 1] + c2
                trace[i, j] = 2

    return backtrace(trace)


def dtw_standard(
    series_a: np.ndarray, series_b: np.ndarray
) -> Tuple[float, np.ndarray]:
    """
    Computes the standard Dynamic Time Warping (DTW) distance and path between two sequences.

    Parameters:
        series_a (np.ndarray): First time series of shape (N,).
        series_b (np.ndarray): Second time series of shape (M,).

    Returns:
        Tuple[float, np.ndarray]: DTW distance and the optimal warping path.
    """
    if series_a.ndim != 1 or series_b.ndim != 1:
        raise ValueError("Input series must be one-dimensional arrays.")

    N, M = len(series_a), len(series_b)
    cost_matrix = np.abs(series_a[:, np.newaxis] - series_b[np.newaxis, :]).astype(
        np.float32
    )

    path = dtw_cpu(cost_matrix)
    distance = cost_matrix[path[0], path[1]].sum()

    return distance, path


def fast_dtw(
    series_a: np.ndarray,
    series_b: np.ndarray,
    radius: int = 1,
    max_dist: Optional[float] = None
) -> Tuple[float, np.ndarray]:
    """
    Computes an approximate DTW distance and path using the FastDTW algorithm.

    Parameters:
        series_a (np.ndarray): First time series.
        series_b (np.ndarray): Second time series.
        radius (int): Radius parameter for FastDTW.
        max_dist (Optional[float]): Maximum distance threshold.

    Returns:
        Tuple[float, np.ndarray]: Approximate DTW distance and warping path.
    """
    def _fast_dtw(a: np.ndarray, b: np.ndarray, radius: int) -> Tuple[float, List[Tuple[int, int]]]:
        if len(a) == 0 or len(b) == 0:
            raise ValueError("Input series must have non-zero length.")

        if len(a) <= radius or len(b) <= radius:
            distance, path = dtw_standard(a, b)
            return distance, list(zip(path[0], path[1]))

        # Downsample
        a_shrink = a[::2]
        b_shrink = b[::2]
        distance, path = _fast_dtw(a_shrink, b_shrink, radius)

        window = expand_window(path, len(a), len(b), radius)
        distance_matrix = np.abs(a[:, np.newaxis] - b[np.newaxis, :]).astype(np.float32)
        distance, path = dtw_cpu_with_window(distance_matrix, window)

        return distance, list(zip(path[0], path[1]))

    def expand_window(
        path: List[Tuple[int, int]], len_a: int, len_b: int, radius: int
    ) -> set:
        window = set()
        for i, j in path:
            for di in range(-radius, radius + 1):
                for dj in range(-radius, radius + 1):
                    ni = 2 * i + di
                    nj = 2 * j + dj
                    if 0 <= ni < len_a and 0 <= nj < len_b:
                        window.add((ni, nj))
        return window

    def dtw_cpu_with_window(x: np.ndarray, window: set) -> Tuple[float, np.ndarray]:
        N, M = x.shape
        cost = np.full((N + 1, M + 1), np.inf, dtype=np.float32)
        trace = -np.ones((N + 1, M + 1), dtype=np.int32)
        cost[0, 0] = 0

        for i in prange(1, N + 1):
            for j in range(1, M + 1):
                if (i - 1, j - 1) in window:
                    c0 = cost[i - 1, j - 1]
                    c1 = cost[i - 1, j]
                    c2 = cost[i, j - 1]

                    if c0 < c1 and c0 < c2:
                        cost[i, j] = x[i - 1, j - 1] + c0
                        trace[i, j] = 0
                    elif c1 < c0 and c1 < c2:
                        cost[i, j] = x[i - 1, j - 1] + c1
                        trace[i, j] = 1
                    else:
                        cost[i, j] = x[i - 1, j - 1] + c2
                        trace[i, j] = 2

        path = backtrace(trace)
        distance = cost[N, M]
        return distance, path

    distance, path = _fast_dtw(series_a, series_b, radius)
    return distance, np.array(path).T


def subsequence_dtw(
    series_a: np.ndarray, series_b: np.ndarray
) -> Tuple[float, np.ndarray]:
    """
    Computes the Subsequence DTW distance and path, finding the best matching subsequence.

    Parameters:
        series_a (np.ndarray): Larger time series.
        series_b (np.ndarray): Subsequence to find within series_a.

    Returns:
        Tuple[float, np.ndarray]: Minimum DTW distance and the warping path.
    """
    if series_a.ndim != 1 or series_b.ndim != 1:
        raise ValueError("Input series must be one-dimensional arrays.")

    N, M = len(series_a), len(series_b)
    if M > N:
        raise ValueError("Subsequence length must be less than or equal to the main series length.")

    min_distance = np.inf
    best_path = np.empty((2, 0), dtype=int)

    for i in range(N - M + 1):
        window_a = series_a[i : i + M]
        distance, path = dtw_standard(window_a, series_b)
        if distance < min_distance:
            min_distance = distance
            best_path = path + i  # Offset path indices by i

    return min_distance, best_path


@njit
def weighted_dtw_cpu(x: np.ndarray, weight: np.ndarray) -> np.ndarray:
    """
    Computes the Weighted DTW alignment path on CPU.

    Parameters:
        x (np.ndarray): Cost matrix.
        weight (np.ndarray): Weight matrix.

    Returns:
        np.ndarray: Optimal warping path.
    """
    N, M = x.shape
    cost = np.full((N + 1, M + 1), np.inf, dtype=np.float32)
    trace = -np.ones((N + 1, M + 1), dtype=np.int32)

    cost[0, 0] = 0

    for i in prange(1, N + 1):
        for j in range(1, M + 1):
            c0 = cost[i - 1, j - 1] + weight[i - 1, j - 1]
            c1 = cost[i - 1, j] + weight[i - 1, j - 1]
            c2 = cost[i, j - 1] + weight[i - 1, j - 1]

            if c0 < c1 and c0 < c2:
                cost[i, j] = x[i - 1, j - 1] + c0
                trace[i, j] = 0
            elif c1 < c0 and c1 < c2:
                cost[i, j] = x[i - 1, j - 1] + c1
                trace[i, j] = 1
            else:
                cost[i, j] = x[i - 1, j - 1] + c2
                trace[i, j] = 2

    return backtrace(trace)


def dtw_weighted(
    series_a: np.ndarray, series_b: np.ndarray, weight_function: Optional[str] = None
) -> Tuple[float, np.ndarray]:
    """
    Computes the Weighted DTW (WDTW) distance and path between two sequences.

    Parameters:
        series_a (np.ndarray): First time series.
        series_b (np.ndarray): Second time series.
        weight_function (Optional[str]): Type of weighting ('linear', 'quadratic', etc.).

    Returns:
        Tuple[float, np.ndarray]: Weighted DTW distance and warping path.
    """
    if series_a.ndim != 1 or series_b.ndim != 1:
        raise ValueError("Input series must be one-dimensional arrays.")

    N, M = len(series_a), len(series_b)
    cost_matrix = np.abs(series_a[:, np.newaxis] - series_b[np.newaxis, :]).astype(
        np.float32
    )

    # Define a simple weighting function (e.g., linear)
    if weight_function == "linear":
        weight = np.arange(N)[:, None] + np.arange(M)[None, :]
    elif weight_function == "quadratic":
        weight = (np.arange(N)[:, None] ** 2) + (np.arange(M)[None, :] ** 2)
    else:
        weight = np.ones_like(cost_matrix, dtype=np.float32)

    path = weighted_dtw_cpu(cost_matrix, weight)
    distance = (cost_matrix[path[0], path[1]] * weight[path[0], path[1]]).sum()

    return distance, path


@njit(parallel=True)
def dtw_constrained_cpu(
    x: np.ndarray, window: int
) -> np.ndarray:
    """
    Computes the Constrained DTW alignment path on CPU using the Sakoe-Chiba band.

    Parameters:
        x (np.ndarray): Cost matrix.
        window (int): Constraint window size.

    Returns:
        np.ndarray: Optimal warping path.
    """
    N, M = x.shape
    cost = np.full((N + 1, M + 1), np.inf, dtype=np.float32)
    trace = -np.ones((N + 1, M + 1), dtype=np.int32)
    cost[0, 0] = 0

    for i in prange(1, N + 1):
        start_j = max(1, i - window)
        end_j = min(M + 1, i + window + 1)
        for j in range(start_j, end_j):
            c0 = cost[i - 1, j - 1]
            c1 = cost[i - 1, j]
            c2 = cost[i, j - 1]

            if c0 < c1 and c0 < c2:
                cost[i, j] = x[i - 1, j - 1] + c0
                trace[i, j] = 0
            elif c1 < c0 and c1 < c2:
                cost[i, j] = x[i - 1, j - 1] + c1
                trace[i, j] = 1
            else:
                cost[i, j] = x[i - 1, j - 1] + c2
                trace[i, j] = 2

    return backtrace(trace)


def dtw_constrained(
    series_a: np.ndarray, series_b: np.ndarray, window: int
) -> Tuple[float, np.ndarray]:
    """
    Computes the Constrained DTW (CDTW) distance and path between two sequences using a Sakoe-Chiba band.

    Parameters:
        series_a (np.ndarray): First time series.
        series_b (np.ndarray): Second time series.
        window (int): Size of the Sakoe-Chiba band.

    Returns:
        Tuple[float, np.ndarray]: Constrained DTW distance and warping path.
    """
    if series_a.ndim != 1 or series_b.ndim != 1:
        raise ValueError("Input series must be one-dimensional arrays.")

    N, M = len(series_a), len(series_b)
    if window < abs(N - M):
        window = abs(N - M)

    cost_matrix = np.abs(series_a[:, np.newaxis] - series_b[np.newaxis, :]).astype(
        np.float32
    )

    path = dtw_constrained_cpu(cost_matrix, window)
    distance = cost_matrix[path[0], path[1]].sum()

    return distance, path

import numpy as np

# Sample time series
series_a = np.array([1, 2, 3, 4, 2, 1])
series_b = np.array([1, 2, 2, 3, 4, 2, 1])

# Standard DTW
distance, path = dtw_standard(series_a, series_b)
print(f"Standard DTW Distance: {distance}")
print(f"Path: {path}")

# FastDTW
distance_fast, path_fast = fast_dtw(series_a, series_b, radius=1)
print(f"FastDTW Distance: {distance_fast}")
print(f"Path: {path_fast}")

# Subsequence DTW
distance_sub, path_sub = subsequence_dtw(series_a, series_b)
print(f"Subsequence DTW Distance: {distance_sub}")
print(f"Path: {path_sub}")

# Weighted DTW with linear weighting
distance_wdtw, path_wdtw = dtw_weighted(series_a, series_b, weight_function='linear')
print(f"Weighted DTW Distance: {distance_wdtw}")
print(f"Path: {path_wdtw}")

# Constrained DTW with Sakoe-Chiba band
distance_cdtw, path_cdtw = dtw_constrained(series_a, series_b, window=1)
print(f"Constrained DTW Distance: {distance_cdtw}")
print(f"Path: {path_cdtw}")