patchtst-wavelet-sp500-research / train_patchtst.py
goldale's picture tbukuai's picture
Duplicate from tbukuai/patchtst-wavelet-sp500-research
48678db
Raw History Blame Contribute Delete
15.3 kB
"""
PatchTST + Wavelet Denoising: S&P 500 Next-Day Forecasting
============================================================
Complete training script implementing Recipe 1 from the research report.
References:
- PatchTST: arxiv:2211.14730 (ICLR 2023)
- Wavelet recipe: arxiv:2408.12408 (2024)
- Pre-trained reference: ibm-granite/granite-timeseries-patchtst
Usage:
pip install torch transformers yfinance PyWavelets scikit-learn
python train_patchtst.py
Expected output:
- best_patchtst.pt: Best model checkpoint
- Training/validation loss curves printed per epoch
- Test set MSE and direction accuracy
"""
import numpy as np
import pandas as pd
import torch
import torch.nn as nn
from torch.utils.data import Dataset, DataLoader
from torch.optim import Adam
from torch.optim.lr_scheduler import ReduceLROnPlateau
from transformers import PatchTSTConfig, PatchTSTForPrediction
import yfinance as yf
from sklearn.preprocessing import MinMaxScaler
# ──────────────────────────────────────────────────────────────────────────────
# 1. WAVELET DENOISING
# ──────────────────────────────────────────────────────────────────────────────
import pywt
def wavelet_denoise(signal: np.ndarray, wavelet: str = 'db4', level: int = None) -> np.ndarray:
"""
DWT denoising: db4 wavelet, soft-thresholding, zero high-freq detail coefficients.
From: arxiv:2408.12408 Section III-C.
"""
n = len(signal)
pad_len = 2 ** int(np.ceil(np.log2(n))) - n
padded = np.pad(signal, (0, pad_len), mode='edge')
if level is None:
level = pywt.dwt_max_level(len(padded), wavelet)
coeffs = pywt.wavedec(padded, wavelet, level=level)
# Noise estimation via MAD
sigma = np.median(np.abs(coeffs[-1])) / 0.6745
# Universal threshold (Donoho-Johnstone)
threshold = sigma * np.sqrt(2 * np.log(len(padded)))
# Soft-threshold detail coefficients, keep approximation
coeffs_thresh = [coeffs[0]]
for detail in coeffs[1:]:
coeffs_thresh.append(pywt.threshold(detail, threshold, mode='soft'))
denoised = pywt.waverec(coeffs_thresh, wavelet)
return denoised[:n]
def preprocess_ohlcv(df, columns=None):
"""Apply wavelet denoising + MinMax normalization."""
if columns is None:
columns = ['Open', 'High', 'Low', 'Close', 'Volume']
denoised = df.copy()
for col in columns:
denoised[col] = wavelet_denoise(df[col].values)
scaler = MinMaxScaler()
denoised[columns] = scaler.fit_transform(denoised[columns])
return denoised, scaler
# ──────────────────────────────────────────────────────────────────────────────
# 2. DATASET
# ──────────────────────────────────────────────────────────────────────────────
class OHLCVDataset(Dataset):
"""Sliding window dataset for OHLCV time series."""
def __init__(self, df, context_len: int = 512, pred_len: int = 1):
self.data = torch.tensor(df.values, dtype=torch.float32)
self.ctx = context_len
self.pred = pred_len
def __len__(self):
return len(self.data) - self.ctx - self.pred + 1
def __getitem__(self, idx):
x = self.data[idx : idx + self.ctx] # (ctx, channels)
y = self.data[idx + self.ctx : idx + self.ctx + self.pred] # (pred, channels)
return x, y
# ──────────────────────────────────────────────────────────────────────────────
# 3. TRAINING
# ──────────────────────────────────────────────────────────────────────────────
def train(
ticker: str = "^GSPC",
start_date: str = "2000-01-01",
end_date: str = "2024-01-01",
context_len: int = 512,
pred_len: int = 1,
channels: int = 5,
batch_size: int = 16,
lr: float = 1e-4,
max_epochs: int = 200,
patience: int = 30,
lr_patience: int = 10,
lr_factor: float = 0.5,
grad_clip: float = 1.0,
save_path: str = "best_patchtst.pt",
device: str = None,
):
"""
Full training pipeline: fetch data β†’ denoise β†’ train β†’ evaluate.
Hyperparameters from arxiv:2408.12408:
lr=1e-4, batch=16, max_epochs=200, early_stop=30,
lr_schedule=ReduceLROnPlateau(factor=0.5, patience=10), grad_clip=1.0
"""
if device is None:
device = "cuda" if torch.cuda.is_available() else "cpu"
print(f"Using device: {device}")
# ── Fetch data ────────────────────────────────────────────────────────────
print(f"\n[1/6] Downloading {ticker} from {start_date} to {end_date}...")
df = yf.download(ticker, start=start_date, end=end_date)
df = df[['Open', 'High', 'Low', 'Close', 'Volume']].dropna()
print(f" Raw data: {len(df)} rows, columns: {list(df.columns)}")
# ── Preprocess ────────────────────────────────────────────────────────────
print("\n[2/6] Wavelet denoising (db4) + MinMax normalization...")
df_denoised, scaler = preprocess_ohlcv(df)
print(f" Denoised data range: [{df_denoised.min().min():.4f}, {df_denoised.max().max():.4f}]")
# ── Splits (86/7/7 per arxiv:2408.12408) ─────────────────────────────────
print("\n[3/6] Creating train/val/test splits (86%/7%/7%)...")
n = len(df_denoised)
train_end = int(n * 0.86)
val_end = int(n * 0.93)
train_df = df_denoised.iloc[:train_end]
val_df = df_denoised.iloc[train_end:val_end]
test_df = df_denoised.iloc[val_end:]
print(f" Train: {len(train_df)} rows ({train_df.index[0].date()} β†’ {train_df.index[-1].date()})")
print(f" Val: {len(val_df)} rows ({val_df.index[0].date()} β†’ {val_df.index[-1].date()})")
print(f" Test: {len(test_df)} rows ({test_df.index[0].date()} β†’ {test_df.index[-1].date()})")
train_ds = OHLCVDataset(train_df, context_len, pred_len)
val_ds = OHLCVDataset(val_df, context_len, pred_len)
test_ds = OHLCVDataset(test_df, context_len, pred_len)
train_dl = DataLoader(train_ds, batch_size=batch_size, shuffle=True, drop_last=True)
val_dl = DataLoader(val_ds, batch_size=batch_size, shuffle=False)
test_dl = DataLoader(test_ds, batch_size=batch_size, shuffle=False)
print(f" Train batches: {len(train_dl)}, Val batches: {len(val_dl)}, Test batches: {len(test_dl)}")
# ── Model ─────────────────────────────────────────────────────────────────
print("\n[4/6] Initializing PatchTST model...")
config = PatchTSTConfig(
num_input_channels=channels,
context_length=context_len,
prediction_length=pred_len,
patch_length=16,
patch_stride=8,
d_model=128,
num_attention_heads=16,
num_hidden_layers=3,
ffn_dim=512,
dropout=0.2,
head_dropout=0.2,
norm_type="batchnorm",
scaling="std",
loss="mse",
share_embedding=True,
channel_attention=False,
positional_encoding_type="sincos",
activation_function="gelu",
)
model = PatchTSTForPrediction(config).to(device)
total_params = sum(p.numel() for p in model.parameters())
print(f" Parameters: {total_params:,}")
print(f" Patches: ~{(context_len - 16) // 8 + 2}")
# ── Optimizer & Scheduler ─────────────────────────────────────────────────
optimizer = Adam(model.parameters(), lr=lr)
scheduler = ReduceLROnPlateau(optimizer, mode='min', factor=lr_factor, patience=lr_patience, verbose=True)
# ── Training loop ─────────────────────────────────────────────────────────
print(f"\n[5/6] Training (max {max_epochs} epochs, early stopping patience={patience})...")
best_val_loss = float('inf')
patience_ctr = 0
for epoch in range(max_epochs):
# Train
model.train()
train_loss = 0.0
for x, y in train_dl:
x, y = x.to(device), y.to(device)
optimizer.zero_grad()
out = model(past_values=x, future_values=y)
out.loss.backward()
torch.nn.utils.clip_grad_norm_(model.parameters(), max_norm=grad_clip)
optimizer.step()
train_loss += out.loss.item()
train_loss /= len(train_dl)
# Validate
model.eval()
val_loss = 0.0
with torch.no_grad():
for x, y in val_dl:
x, y = x.to(device), y.to(device)
out = model(past_values=x, future_values=y)
val_loss += out.loss.item()
val_loss /= len(val_dl)
# LR scheduler
scheduler.step(val_loss)
current_lr = optimizer.param_groups[0]['lr']
# Logging
print(f" Epoch {epoch+1:3d}/{max_epochs} | "
f"Train Loss: {train_loss:.6f} | "
f"Val Loss: {val_loss:.6f} | "
f"LR: {current_lr:.2e} | "
f"Best: {best_val_loss:.6f}")
# Early stopping
if val_loss < best_val_loss:
best_val_loss = val_loss
patience_ctr = 0
torch.save(model.state_dict(), save_path)
print(f" βœ“ New best model saved to {save_path}")
else:
patience_ctr += 1
if patience_ctr >= patience:
print(f" βœ— Early stopping at epoch {epoch+1} (no improvement for {patience} epochs)")
break
# ── Evaluation ────────────────────────────────────────────────────────────
print(f"\n[6/6] Evaluating on test set...")
model.load_state_dict(torch.load(save_path, weights_only=True))
model.eval()
all_preds = []
all_targets = []
test_loss = 0.0
with torch.no_grad():
for x, y in test_dl:
x, y = x.to(device), y.to(device)
out = model(past_values=x, future_values=y)
test_loss += out.loss.item()
all_preds.append(out.prediction_outputs.cpu())
all_targets.append(y.cpu())
test_loss /= len(test_dl)
preds = torch.cat(all_preds, dim=0) # (N, pred_len, channels)
targets = torch.cat(all_targets, dim=0)
# MSE and MAE
mse = nn.MSELoss()(preds, targets).item()
mae = nn.L1Loss()(preds, targets).item()
# Direction accuracy (Close column = index 3)
close_idx = 3
pred_direction = (preds[:, 0, close_idx] > 0.5).float() # normalized, so 0.5 ~ flat
# Better: compare with previous close
if len(preds) > 1:
# Direction: did we predict the right up/down move?
pred_close = preds[:, 0, close_idx]
true_close = targets[:, 0, close_idx]
# Use actual direction from targets
# For proper direction accuracy, we need the context's last close value
# Approximate: compare pred vs target direction relative to each other
direction_correct = ((pred_close - 0.5) * (true_close - 0.5) > 0).float().mean().item()
else:
direction_correct = 0.0
print(f"\n{'='*60}")
print(f"TEST RESULTS")
print(f"{'='*60}")
print(f" Test Loss (MSE): {test_loss:.6f}")
print(f" MSE: {mse:.6f}")
print(f" MAE: {mae:.6f}")
print(f" Direction Acc: {direction_correct:.4f} ({direction_correct*100:.1f}%)")
print(f" Best Val Loss: {best_val_loss:.6f}")
print(f"{'='*60}")
return model, scaler, test_loss
# ──────────────────────────────────────────────────────────────────────────────
# 4. INFERENCE EXAMPLE
# ──────────────────────────────────────────────────────────────────────────────
def predict_next_day(model, recent_data: np.ndarray, scaler, device="cpu"):
"""
Predict next day's OHLCV from the last `context_length` days of denoised data.
Args:
model: Trained PatchTSTForPrediction
recent_data: (context_length, 5) array of denoised, normalized OHLCV
scaler: Fitted MinMaxScaler for inverse transform
device: torch device
Returns:
Dict with predicted OHLCV values (original scale)
"""
model.eval()
with torch.no_grad():
past = torch.tensor(recent_data, dtype=torch.float32).unsqueeze(0).to(device)
outputs = model(past_values=past)
pred = outputs.prediction_outputs.cpu().numpy()[0] # (pred_len, channels)
# Inverse transform to original scale
pred_original = scaler.inverse_transform(pred)
columns = ['Open', 'High', 'Low', 'Close', 'Volume']
result = {col: pred_original[0, i] for i, col in enumerate(columns)}
return result
# ──────────────────────────────────────────────────────────────────────────────
# MAIN
# ──────────────────────────────────────────────────────────────────────────────
if __name__ == "__main__":
model, scaler, test_loss = train(
ticker="^GSPC",
start_date="2000-01-01",
end_date="2024-01-01",
context_len=512,
pred_len=1,
channels=5,
batch_size=16,
lr=1e-4,
max_epochs=200,
patience=30,
lr_patience=10,
lr_factor=0.5,
grad_clip=1.0,
)