"""Price and Volume Feature Families — core feature families derived from a single asset's price and volume history."""

from __future__ import annotations

import warnings
from datetime import datetime

import numpy as np
import plotly.graph_objects as go
import polars as pl
from plotly.subplots import make_subplots

from utils.paths import get_chapter_dir
from utils.style import COLORS  # importing utils.style sets the ml4t Plotly template as default

warnings.filterwarnings("ignore")

SEED = 42
START_DATE = "2015-01-01"

from ml4t.engineer import compute_features
from ml4t.engineer.core.registry import get_registry

registry = get_registry()
all_features = registry.list_all()

# Features by category
categories = {}
for name in all_features:
    metadata = registry.get(name)
    categories.setdefault(metadata.category, []).append(name)

print(f"Total features available: {len(all_features)}\n")
print("Features by category:")
for cat, feats in sorted(categories.items()):
    examples = ", ".join(feats[:4])
    suffix = ", ..." if len(feats) > 4 else ""
    print(f"  {cat:20s}: {len(feats):3d}  ({examples}{suffix})")

rsi_meta = registry.get("rsi")
print(f"Feature:     {rsi_meta.name}")
print(f"Category:    {rsi_meta.category}")
print(f"Formula:     {rsi_meta.formula}")
print(f"Parameters:  {rsi_meta.parameters}")
print(f"Input type:  {rsi_meta.input_type}")

from data import load_etfs

spy_quick = load_etfs().filter(pl.col("symbol") == "SPY").sort("timestamp").tail(500)

result = compute_features(spy_quick, ["rsi", "atr", "sma"])
new_cols = [c for c in result.columns if c not in spy_quick.columns]
print(f"Computed {len(new_cols)} feature columns: {new_cols}")
result.select(["timestamp", "close"] + new_cols).tail(5)

# Load and sort by symbol, then timestamp (CRITICAL for .over() operations)
etfs = load_etfs().sort(["symbol", "timestamp"])

# Filter date range
etfs = etfs.filter(pl.col("timestamp") >= datetime.strptime(START_DATE, "%Y-%m-%d"))

# For single-asset demos: SPY
spy = etfs.filter(pl.col("symbol") == "SPY").sort("timestamp")

# For cross-sectional demos: subset of liquid ETFs
cs_symbols = ["SPY", "QQQ", "IWM", "TLT", "GLD", "XLF", "XLE", "XLK"]
cs_etfs = etfs.filter(pl.col("symbol").is_in(cs_symbols)).sort(["symbol", "timestamp"])

print(f"SPY: {len(spy):,} rows")
print(f"Cross-sectional universe: {len(cs_etfs):,} rows, {cs_etfs['symbol'].n_unique()} symbols")
print(f"Date range: {spy['timestamp'].min()} to {spy['timestamp'].max()}")

def compute_returns_manual(df: pl.DataFrame, horizons: list[int]) -> pl.DataFrame:
    """
    Compute returns for multiple horizons using pure Polars.

    Note: This is for teaching. Use ml4t-engineer in production.
    """
    return_exprs = []

    for h in horizons:
        # Simple returns
        return_exprs.append(pl.col("close").pct_change(h).alias(f"ret_{h}d"))
        # Log returns
        return_exprs.append(
            (pl.col("close").log() - pl.col("close").log().shift(h)).alias(f"logret_{h}d")
        )

    return df.with_columns(return_exprs)


# Single with_columns for efficiency
returns_df = compute_returns_manual(spy, horizons=[1, 5, 21])

print("Return features computed:")
returns_df.select(["timestamp", "close", "ret_1d", "ret_5d", "ret_21d"]).tail(10)

# Skip-1 momentum: shift(1) before computing the horizon return
skip1_df = spy.with_columns(
    [
        # Standard 21-day momentum
        pl.col("close").pct_change(21).alias("mom_21d"),
        # Skip-1: use yesterday's close as numerator
        ((pl.col("close").shift(1) / pl.col("close").shift(21)) - 1).alias("mom_21d_skip1"),
    ]
)

# Show correlation - should be high but not identical
correlation = skip1_df.drop_nulls().select(
    [pl.corr("mom_21d", "mom_21d_skip1").alias("correlation")]
)
print(f"Correlation between standard and skip-1 momentum: {correlation[0, 0]:.4f}")

session_df = spy.with_columns(
    [
        # Overnight (gap)
        ((pl.col("open") / pl.col("close").shift(1)) - 1).alias("overnight_ret"),
        # Intraday
        ((pl.col("close") / pl.col("open")) - 1).alias("intraday_ret"),
        # Total for reference
        pl.col("close").pct_change().alias("total_ret"),
    ]
)

# Verify decomposition: (1 + overnight) * (1 + intraday) ≈ (1 + total)
session_df = session_df.with_columns(
    ((1 + pl.col("overnight_ret")) * (1 + pl.col("intraday_ret")) - 1).alias("reconstructed")
)

print("Session return decomposition:")
session_df.select(
    ["timestamp", "overnight_ret", "intraday_ret", "total_ret", "reconstructed"]
).tail(10)

from ml4t.engineer.features.volatility import atr

# MA distance scaled by ATR
ma_df = spy.with_columns(
    [
        pl.col("close").rolling_mean(21).alias("sma_21"),
        atr("high", "low", "close", period=21).alias("atr_21"),
    ]
).with_columns(
    [
        # Raw distance (not comparable across assets)
        (pl.col("close") - pl.col("sma_21")).alias("ma_dist_raw"),
        # Vol-scaled (standardized)
        ((pl.col("close") - pl.col("sma_21")) / pl.col("atr_21")).alias("ma_dist_scaled"),
    ]
)

# Visualize MA distance
fig = make_subplots(
    rows=3,
    cols=1,
    shared_xaxes=True,
    subplot_titles=["SPY Close with SMA(21)", "Raw MA Distance", "Vol-Scaled MA Distance"],
    vertical_spacing=0.08,
)

n = 252  # Last year
fig.add_trace(
    go.Scatter(x=ma_df["timestamp"].to_list()[-n:], y=ma_df["close"].to_list()[-n:], name="Close"),
    row=1,
    col=1,
)
fig.add_trace(
    go.Scatter(
        x=ma_df["timestamp"].to_list()[-n:],
        y=ma_df["sma_21"].to_list()[-n:],
        name="SMA(21)",
        line=dict(dash="dash"),
    ),
    row=1,
    col=1,
)
fig.add_trace(
    go.Scatter(
        x=ma_df["timestamp"].to_list()[-n:],
        y=ma_df["ma_dist_raw"].to_list()[-n:],
        name="Raw",
        fill="tozeroy",
    ),
    row=2,
    col=1,
)
fig.add_trace(
    go.Scatter(
        x=ma_df["timestamp"].to_list()[-n:],
        y=ma_df["ma_dist_scaled"].to_list()[-n:],
        name="Vol-Scaled",
        fill="tozeroy",
    ),
    row=3,
    col=1,
)
fig.add_hline(y=0, line_dash="dash", line_color=COLORS["neutral"], row=2, col=1)
fig.add_hline(y=0, line_dash="dash", line_color=COLORS["neutral"], row=3, col=1)

fig.update_yaxes(title_text="Price ($)", row=1, col=1)
fig.update_yaxes(title_text="Distance ($)", row=2, col=1)
fig.update_yaxes(title_text="Distance (ATR units)", row=3, col=1)
fig.update_xaxes(title_text="Date", row=3, col=1)
fig.update_layout(
    height=600, title="Volatility-scaling standardizes MA distance for cross-asset comparison"
)
fig.show()

def rolling_regression_slope(
    df: pl.DataFrame,
    price_col: str = "close",
    period: int = 21,
) -> pl.DataFrame:
    """Rolling OLS slope, normalized by mean price."""
    return df.with_columns(
        pl.col(price_col)
        .rolling_map(
            lambda s: np.polyfit(np.arange(len(s)), s.to_numpy(), 1)[0] / s.mean() * len(s),
            window_size=period,
            min_periods=period,
        )
        .alias(f"slope_{period}d")
    )


slope_df = rolling_regression_slope(spy, period=21)
print("Rolling regression slope (normalized %-change over 21d):")
slope_df.select(["timestamp", "close", "slope_21d"]).tail(10)

reversal_df = spy.with_columns(
    [
        pl.col("close").pct_change(1).alias("ret_1d"),
        (-pl.col("close").pct_change(1)).alias("reversal_1d"),
        # Z-scored reversal for comparability
        (
            (-pl.col("close").pct_change(1) - (-pl.col("close").pct_change(1)).rolling_mean(63))
            / (-pl.col("close").pct_change(1)).rolling_std(63)
        ).alias("reversal_zscore"),
    ]
)

print("Short-term reversal features:")
reversal_df.select(["timestamp", "close", "ret_1d", "reversal_1d"]).tail(10)

dist_ma_df = spy.with_columns(
    [
        pl.col("close").rolling_mean(50).alias("ma_50"),
        atr("high", "low", "close", period=14).alias("atr_14"),
    ]
).with_columns(
    ((pl.col("close") - pl.col("ma_50")) / pl.col("atr_14").clip(1e-10, None)).alias(
        "dist_to_ma_atr"
    ),
)

print("Distance to MA in ATR units:")
dist_ma_df.select(["timestamp", "close", "ma_50", "atr_14", "dist_to_ma_atr"]).tail(10)

from ml4t.engineer.features.volatility import realized_volatility, yang_zhang_volatility

vol_df = spy.with_columns(
    pl.col("close").pct_change().alias("ret"),
).with_columns(
    [
        # Close-to-close (annualized) — note: realized_volatility expects returns, not prices
        realized_volatility("ret", period=21, annualize=True).alias("vol_cc_21"),
        # Yang-Zhang (most efficient)
        yang_zhang_volatility("open", "high", "low", "close", period=21, annualize=True).alias(
            "vol_yz_21"
        ),
    ]
)

print("Volatility comparison:")
vol_df.select(["timestamp", "close", "vol_cc_21", "vol_yz_21"]).tail(10)

# Manual implementations for teaching
def parkinson_vol(period: int = 21) -> pl.Expr:
    """Parkinson range-based volatility (annualized)."""
    log_hl_sq = ((pl.col("high").log() - pl.col("low").log()) ** 2) / (4 * np.log(2))
    return (log_hl_sq.rolling_mean(period) * 252).sqrt()

def garman_klass_vol(period: int = 21) -> pl.Expr:
    """Garman-Klass OHLC volatility (annualized)."""
    log_hl_sq = 0.5 * (pl.col("high").log() - pl.col("low").log()) ** 2
    log_co_sq = (2 * np.log(2) - 1) * (pl.col("close").log() - pl.col("open").log()) ** 2
    return ((log_hl_sq - log_co_sq).rolling_mean(period) * 252).sqrt()

# Compare all four estimators
vol_compare_df = spy.with_columns(
    pl.col("close").pct_change().alias("ret"),
).with_columns(
    [
        # Close-to-close — realized_volatility expects returns, not prices
        realized_volatility("ret", period=21, annualize=True).alias("vol_cc"),
        # Parkinson (manual — works on OHLC directly)
        parkinson_vol(period=21).alias("vol_parkinson"),
        # Garman-Klass (manual — works on OHLC directly)
        garman_klass_vol(period=21).alias("vol_gk"),
        # Yang-Zhang (library — works on OHLC directly, most efficient)
        yang_zhang_volatility("open", "high", "low", "close", period=21, annualize=True).alias(
            "vol_yz"
        ),
    ]
)

# 4-panel comparison
n = 504  # Last ~2 years for visualization
fig = make_subplots(
    rows=2,
    cols=2,
    shared_xaxes=True,
    shared_yaxes=True,
    subplot_titles=["Close-to-Close", "Parkinson (H-L)", "Garman-Klass (OHLC)", "Yang-Zhang"],
    vertical_spacing=0.10,
    horizontal_spacing=0.06,
)

for idx, (col, name) in enumerate(
    [("vol_cc", "CC"), ("vol_parkinson", "Park"), ("vol_gk", "GK"), ("vol_yz", "YZ")]
):
    row, c = idx // 2 + 1, idx % 2 + 1
    fig.add_trace(
        go.Scatter(
            x=vol_compare_df["timestamp"].to_list()[-n:],
            y=vol_compare_df[col].to_list()[-n:],
            name=name,
        ),
        row=row,
        col=c,
    )

fig.update_yaxes(title_text="Annualized volatility", col=1)
fig.update_xaxes(title_text="Date", row=2)
fig.update_layout(
    height=500,
    title="Range-based estimators track SPY volatility more smoothly than close-to-close (21-day)",
)
fig.show()

import matplotlib.pyplot as plt

# Prepare data — use full history for a window that includes a volatility spike
vol_plot = vol_compare_df.drop_nulls(subset=["vol_cc", "vol_parkinson", "vol_gk", "vol_yz"])

fig_mpl, ax = plt.subplots(figsize=(12, 5))

# SPY close price on secondary axis (context for vol spikes)
ax2 = ax.twinx()
dates = vol_plot["timestamp"].to_list()
ax2.fill_between(dates, vol_plot["close"].to_list(), alpha=0.08, color=COLORS["neutral"])
ax2.set_ylabel("SPY Close ($)", color=COLORS["neutral"])
ax2.tick_params(axis="y", labelcolor=COLORS["neutral"])

# Vol estimators on primary axis (convert to percentage for readability).
# Line styles vary alongside color so the four series stay distinct in grayscale print.
estimators = [
    ("vol_cc", "Close-to-Close (eff. 1.0\u00d7)", "-", COLORS["blue"]),
    ("vol_parkinson", "Parkinson (eff. 5.2\u00d7)", "--", COLORS["amber"]),
    ("vol_gk", "Garman-Klass (eff. 7.4\u00d7)", ":", COLORS["copper"]),
    ("vol_yz", "Yang-Zhang (eff. 8.4\u00d7)", "-.", COLORS["slate"]),
]

for col, label, ls, color in estimators:
    vals = [v * 100 for v in vol_plot[col].to_list()]
    ax.plot(dates, vals, ls=ls, color=color, label=label, linewidth=1.2)

ax.set_ylabel("Annualized Volatility (%)")
ax.set_title("Rolling Volatility Estimates: Four OHLC Estimators on SPY")
ax.legend(frameon=False, fontsize=8, loc="upper left")
ax.set_zorder(ax2.get_zorder() + 1)
ax.patch.set_visible(False)

plt.show()

# Persist source data so the book figure script can re-render at print resolution
# without re-executing this notebook.
_FIG_8_3_ARTIFACT = (
    get_chapter_dir(8)
    / "output"
    / "book_figure_artifacts"
    / "figure_8_3_ohlc_volatility_estimators.parquet"
)
_FIG_8_3_ARTIFACT.parent.mkdir(parents=True, exist_ok=True)
vol_plot.select(
    ["timestamp", "close", "vol_cc", "vol_parkinson", "vol_gk", "vol_yz"]
).write_parquet(_FIG_8_3_ARTIFACT)

from ml4t.engineer.features.volatility import volatility_of_volatility

vov_df = spy.with_columns(
    [
        yang_zhang_volatility("open", "high", "low", "close", period=21, annualize=True).alias(
            "vol"
        ),
        volatility_of_volatility("close", vol_period=21, vov_period=21).alias("vol_of_vol"),
    ]
)

print("Vol-of-vol (last 10 rows):")
vov_df.select(["timestamp", "vol", "vol_of_vol"]).tail(10)

SHORT_VOL_WINDOW = 10
LONG_VOL_WINDOW = 63
PERCENTILE_LOOKBACK = 252

vol_state_df = (
    spy.with_columns(
        pl.col("close").log().diff().alias("log_return"),
    )
    .with_columns(
        [
            (pl.col("log_return").rolling_std(SHORT_VOL_WINDOW) * np.sqrt(252)).alias("vol_short"),
            (pl.col("log_return").rolling_std(LONG_VOL_WINDOW) * np.sqrt(252)).alias("vol_long"),
        ]
    )
    .with_columns(
        # Vol ratio: short/long — >1 means expansion, <1 contraction
        (pl.col("vol_short") / pl.col("vol_long").clip(1e-10, None)).alias("vol_ratio"),
    )
)

# Vol percentile: rolling rank over 252 days
vol_state_df = vol_state_df.with_columns(
    pl.col("vol_short")
    .rolling_map(
        lambda s: 100 * (s.rank().last() - 1) / max(len(s) - 1, 1),
        window_size=PERCENTILE_LOOKBACK,
        min_periods=PERCENTILE_LOOKBACK // 2,
    )
    .alias("vol_percentile")
)

# Vol decile: clipped 0-9
vol_state_df = vol_state_df.with_columns(
    (pl.col("vol_percentile").fill_nan(None).fill_null(-10.0) / 10)
    .floor()
    .clip(0, 9)
    .cast(pl.Int32)
    .alias("vol_decile")
)

print("Volatility state features:")
vol_state_df.select(["timestamp", "vol_ratio", "vol_percentile", "vol_decile"]).tail(10)

from ml4t.engineer.features.regime import fractal_efficiency, trend_intensity_index, variance_ratio

vr_exprs = variance_ratio("close", periods=[5], window=21)
regime_ind_df = spy.with_columns(
    [
        vr_exprs["vr_5"].alias("variance_ratio_21d"),
        fractal_efficiency("close", period=21).alias("fractal_eff_21d"),
        trend_intensity_index("close", period=21).alias("trend_intensity_21d"),
    ]
)

print("Regime indicators:")
regime_ind_df.select(
    ["timestamp", "variance_ratio_21d", "fractal_eff_21d", "trend_intensity_21d"]
).tail(10)

volume_df = spy.with_columns(
    [
        (pl.col("volume") * pl.col("close")).alias("dollar_volume"),
        # Log transform for better distribution
        (pl.col("volume") * pl.col("close")).log().alias("log_dollar_volume"),
    ]
)

print("Dollar volume:")
volume_df.select(["timestamp", "close", "volume", "dollar_volume"]).tail(10)

rel_vol_df = spy.with_columns(
    [
        (pl.col("volume") / pl.col("volume").rolling_mean(21)).alias("rel_volume"),
        # Volume shock: z-score of log volume
        (
            (pl.col("volume").log() - pl.col("volume").log().rolling_mean(21))
            / pl.col("volume").log().rolling_std(21)
        ).alias("volume_zscore"),
    ]
)

# Visualize relative volume features
fig = make_subplots(
    rows=3,
    cols=1,
    shared_xaxes=True,
    subplot_titles=["SPY Close", "Relative Volume", "Volume Z-Score"],
    vertical_spacing=0.08,
)

fig.add_trace(
    go.Scatter(
        x=rel_vol_df["timestamp"].to_list()[-n:], y=rel_vol_df["close"].to_list()[-n:], name="Close"
    ),
    row=1,
    col=1,
)
fig.add_trace(
    go.Scatter(
        x=rel_vol_df["timestamp"].to_list()[-n:],
        y=rel_vol_df["rel_volume"].to_list()[-n:],
        name="Rel Vol",
        fill="tozeroy",
    ),
    row=2,
    col=1,
)
fig.add_hline(y=1, line_dash="dash", line_color=COLORS["neutral"], row=2, col=1)
fig.add_trace(
    go.Scatter(
        x=rel_vol_df["timestamp"].to_list()[-n:],
        y=rel_vol_df["volume_zscore"].to_list()[-n:],
        name="Vol Z-Score",
    ),
    row=3,
    col=1,
)
fig.add_hline(y=0, line_dash="dash", line_color=COLORS["neutral"], row=3, col=1)
fig.add_hline(y=2, line_dash="dash", line_color=COLORS["negative"], row=3, col=1)
fig.add_hline(y=-2, line_dash="dash", line_color=COLORS["negative"], row=3, col=1)

fig.update_yaxes(title_text="Price ($)", row=1, col=1)
fig.update_yaxes(title_text="Rel. volume (×avg)", row=2, col=1)
fig.update_yaxes(title_text="Volume z-score", row=3, col=1)
fig.update_xaxes(title_text="Date", row=3, col=1)
fig.update_layout(height=600, title="Relative volume and z-scores isolate trading-activity spikes")
fig.show()

# Daily VWAP using typical price as proxy
vwap_df = (
    spy.with_columns(
        [
            # Typical price as VWAP proxy for daily data
            ((pl.col("high") + pl.col("low") + pl.col("close")) / 3).alias("typical_price"),
        ]
    )
    .with_columns(
        [
            # Rolling VWAP (volume-weighted rolling mean)
            (
                (pl.col("typical_price") * pl.col("volume")).rolling_sum(5)
                / pl.col("volume").rolling_sum(5)
            ).alias("vwap_5d"),
        ]
    )
    .with_columns(
        [
            ((pl.col("close") / pl.col("vwap_5d")) - 1).alias("vwap_distance"),
        ]
    )
)

print("VWAP distance:")
vwap_df.select(["timestamp", "close", "vwap_5d", "vwap_distance"]).tail(10)

# Add momentum to cross-sectional data
cs_mom = cs_etfs.with_columns(
    [
        pl.col("close").pct_change(21).over("symbol").alias("mom_21d"),
    ]
)

# Cross-sectional rank within each day
cs_ranked = cs_mom.with_columns(
    [
        # Rank: 0 = lowest momentum, 1 = highest
        (
            pl.col("mom_21d").rank().over("timestamp") / pl.col("symbol").count().over("timestamp")
        ).alias("mom_rank"),
    ]
)

# Show one day
sample_date = cs_ranked["timestamp"].max()
print(f"Cross-sectional ranks for {sample_date}:")
(
    cs_ranked.filter(pl.col("timestamp") == sample_date)
    .select(["symbol", "mom_21d", "mom_rank"])
    .sort("mom_rank", descending=True)
)

# Vol-scaled momentum: return / realized vol
cs_vol_scaled = cs_etfs.with_columns(
    [
        pl.col("close").pct_change(21).over("symbol").alias("mom_21d"),
        (
            pl.col("close").pct_change().over("symbol").rolling_std(21).over("symbol")
            * np.sqrt(252)
        ).alias("vol_21d"),
    ]
).with_columns(
    (pl.col("mom_21d") / pl.col("vol_21d").clip(1e-10, None)).alias("vol_scaled_mom"),
)

# Cross-sectional percentile rank
cs_vol_ranked = cs_vol_scaled.with_columns(
    (
        pl.col("vol_scaled_mom").rank().over("timestamp")
        / pl.col("symbol").count().over("timestamp")
    ).alias("vol_scaled_rank"),
)

# Compare raw vs vol-scaled rankings
print(f"Vol-scaled cross-sectional momentum ({sample_date}):")
(
    cs_vol_ranked.filter(pl.col("timestamp") == sample_date)
    .select(["symbol", "mom_21d", "vol_21d", "vol_scaled_mom", "vol_scaled_rank"])
    .sort("vol_scaled_rank", descending=True)
)

# Z-score within each day
cs_zscored = cs_mom.with_columns(
    [
        # Standard z-score
        (
            (pl.col("mom_21d") - pl.col("mom_21d").mean().over("timestamp"))
            / pl.col("mom_21d").std().over("timestamp")
        ).alias("mom_zscore"),
        # Robust z-score (using median and MAD)
        # MAD × 1.4826 ≈ σ for normal data (1.4826 = 1/Φ⁻¹(3/4))
        (
            (pl.col("mom_21d") - pl.col("mom_21d").median().over("timestamp"))
            / (
                (pl.col("mom_21d") - pl.col("mom_21d").median().over("timestamp"))
                .abs()
                .median()
                .over("timestamp")
                * 1.4826
            )
        ).alias("mom_zscore_robust"),
    ]
)

# Show one day
print(f"Cross-sectional z-scores for {sample_date}:")
(
    cs_zscored.filter(pl.col("timestamp") == sample_date)
    .select(["symbol", "mom_21d", "mom_zscore", "mom_zscore_robust"])
    .sort("mom_zscore", descending=True)
)

from ml4t.engineer.features.risk import (
    conditional_value_at_risk,
    downside_deviation,
    tail_ratio,
    value_at_risk,
)

risk_df = spy.with_columns(
    [
        pl.col("close").pct_change().alias("ret"),
    ]
).with_columns(
    [
        value_at_risk("ret", confidence_level=0.95, window=63).alias("var_5pct_63d"),
        conditional_value_at_risk("ret", confidence_level=0.95, window=63).alias("cvar_5pct_63d"),
        downside_deviation("ret", window=63).alias("downside_dev_63d"),
        tail_ratio("ret", confidence_level=0.95, window=63).alias("tail_ratio_63d"),
    ]
)

print("Risk features (last 10 rows):")
risk_df.select(
    ["timestamp", "var_5pct_63d", "cvar_5pct_63d", "downside_dev_63d", "tail_ratio_63d"]
).tail(10)

from ml4t.engineer.features.fdiff import ffdiff

# Fractional differencing preserves memory while achieving stationarity
# d=0.5 is a common starting point; use find_optimal_d() for data-driven choice
ffd_df = spy.with_columns(
    [
        ffdiff("close", d=0.5).alias("close_ffd_05"),
        ffdiff("close", d=1.0).alias("close_ffd_10"),  # Equivalent to first difference
    ]
)

print("Fractional differencing (d=0.5 vs d=1.0):")
ffd_df.select(["timestamp", "close", "close_ffd_05", "close_ffd_10"]).tail(10)