# SPDX-License-Identifier: MIT """YS-EN-0003 ยท data and image encodings. Data: the UCR/UEA StarLightCurves benchmark (3 classes, 1,000 train / 8,236 test, length 1,024), loaded with aeon. Resample 0 is the official split; resample k>0 is aeon's stratified resample with random_state=k, the protocol behind the published results we compare against. Encodings (all computed on the GPU from the raw series, values kept as floats; no colour map): gasf Gramian Angular Summation Field cos(phi_i + phi_j) gadf Gramian Angular Difference Field sin(phi_i - phi_j) mtf Markov Transition Field (Q quantile bins) rp unthresholded recurrence plot |x_i - x_j| line the curve drawn as a line plot stack [gasf, gadf, mtf] as the three channels """ from __future__ import annotations import os from pathlib import Path import numpy as np try: # encodings need torch; loading data does not import torch except ImportError: # pragma: no cover torch = None WORK = Path(os.environ.get("YS_WORK", Path.home() / "ys-work")) # downloads and caches live here DATA_DIR = WORK / "slcdata" CLASS_NAMES = {"1": "Cepheid", "2": "Eclipsing binary", "3": "RR Lyrae"} def load(resample: int = 0): from aeon.benchmarking.resampling import stratified_resample_data_indices from aeon.datasets import load_classification X, y = load_classification("StarLightCurves", split="train", extract_path=str(DATA_DIR)) Xt, yt = load_classification("StarLightCurves", split="test", extract_path=str(DATA_DIR)) X, Xt = X[:, 0, :].astype(np.float32), Xt[:, 0, :].astype(np.float32) if resample > 0: tr, te = stratified_resample_data_indices(y, yt, random_state=resample) XA, yA = np.concatenate([X, Xt]), np.concatenate([y, yt]) X, y, Xt, yt = XA[tr], yA[tr], XA[te], yA[te] return X, y.astype(int) - 1, Xt, yt.astype(int) - 1 OGLE_TRAIN = {0: 152, 1: 275, 2: 573} # StarLightCurves' training mix def load_ogle(resample: int = 0, with_period: bool = False, variant: str = "", with_ids: bool = False): """OGLE-SLC (build_ogle.py): same size and class mix as StarLightCurves, but every curve keeps its catalogue period. Resample 0 is the saved split; resample k>0 draws a fresh stratified 152 / 275 / 573 training set.""" d = np.load(WORK / "ogle" / f"ogle_slc{variant}.npz") X, y, P, train = d["X"].astype(np.float32), d["y"].astype(int), d["period"], d["train"] if resample > 0: rng = np.random.default_rng(1000 + resample) train = np.zeros(len(y), bool) for c, k in OGLE_TRAIN.items(): train[rng.choice(np.where(y == c)[0], k, replace=False)] = True out = (X[train], y[train], X[~train], y[~train]) if with_period: out += (P[train], P[~train]) return out + (d["ids"][~train],) if with_ids else out def subsample(y: np.ndarray, n: int, seed: int) -> np.ndarray: """Stratified subset of n training indices (for learning curves).""" rng = np.random.default_rng(seed) idx = [] for c in np.unique(y): ci = np.where(y == c)[0] k = max(1, round(n * len(ci) / len(y))) idx.append(rng.choice(ci, k, replace=False)) return np.sort(np.concatenate(idx)) # ---------------------------------------------------------------- encodings def paa(x: torch.Tensor, m: int) -> torch.Tensor: """Piecewise aggregate approximation: (B, L) -> (B, m) segment means.""" return torch.nn.functional.adaptive_avg_pool1d(x[:, None, :], m)[:, 0, :] def minmax(x: torch.Tensor, lo=-1.0, hi=1.0) -> torch.Tensor: a = x.amin(1, keepdim=True) b = x.amax(1, keepdim=True) return lo + (hi - lo) * (x - a) / (b - a).clamp_min(1e-8) def gasf(x, m): z = minmax(paa(x, m)).clamp(-1, 1) s = torch.sqrt((1 - z * z).clamp_min(0)) return z[:, :, None] * z[:, None, :] - s[:, :, None] * s[:, None, :] def gadf(x, m): z = minmax(paa(x, m)).clamp(-1, 1) s = torch.sqrt((1 - z * z).clamp_min(0)) return s[:, :, None] * z[:, None, :] - z[:, :, None] * s[:, None, :] def mtf(x, m, q: int = 8): """Markov transition field: bin the full-resolution series into q quantile bins, estimate the bin-to-bin transition matrix W, set M_ij = W[b_i, b_j], then average-pool the L x L field to m x m.""" B, L = x.shape qs = torch.quantile(x, torch.linspace(0, 1, q + 1, device=x.device)[1:-1], dim=1).T # (B, q-1) b = torch.searchsorted(qs.contiguous(), x.contiguous()) # (B, L) bins 0..q-1 flat = torch.zeros(B, q * q, device=x.device) flat.scatter_add_(1, b[:, :-1] * q + b[:, 1:], torch.ones(B, L - 1, device=x.device)) W = flat.view(B, q, q) W = W / W.sum(2, keepdim=True).clamp_min(1) # Block-averaging the L x L field to m x m equals mixing W by each segment's # bin fractions: M_ij = sum_qr f_i(q) W[q, r] f_j(r). Exact, and avoids L x L. seg = torch.nn.functional.one_hot(b, q).float() # (B, L, q) segp = torch.nn.functional.adaptive_avg_pool1d(seg.transpose(1, 2), m) # (B, q, m) bin weights per segment M = torch.einsum("bqi,bqr,brj->bij", segp, W, segp) return minmax(M.flatten(1)).view(B, m, m) def rp(x, m): z = minmax(paa(x, m), 0, 1) d = (z[:, :, None] - z[:, None, :]).abs() return 1 - 2 * d # 1 on the diagonal, -1 far apart def line(x, m, width: float = 1.5): """Antialiased line plot rendered directly as an m x m intensity image.""" z = minmax(paa(x, m * 2), 0, 1) # 2 samples per column for smoothness rows = torch.arange(m, device=x.device).float() ypix = (1 - z) * (m - 1) # (B, 2m) ypix = torch.nn.functional.avg_pool1d(ypix[:, None, :], 2)[:, 0, :] # (B, m) lo = torch.minimum(ypix, torch.roll(ypix, 1, 1)) hi = torch.maximum(ypix, torch.roll(ypix, 1, 1)) lo[:, 0], hi[:, 0] = ypix[:, 0], ypix[:, 0] r = rows[None, :, None] dist = torch.clamp(lo[:, None, :] - r, min=0) + torch.clamp(r - hi[:, None, :], min=0) img = torch.clamp(1 - dist / width, 0, 1) return img * 2 - 1 ENCODERS = {"gasf": gasf, "gadf": gadf, "mtf": mtf, "rp": rp, "line": line} def encode(x: torch.Tensor, name: str, m: int = 224) -> torch.Tensor: """(B, L) series -> (B, 3, m, m) image in [-1, 1].""" if name == "stack": return torch.stack([gasf(x, m), gadf(x, m), mtf(x, m)], 1) img = ENCODERS[name](x, m) return img[:, None].expand(-1, 3, -1, -1)