"""Preprocessing and normalisation.

Pipeline order (see RESEARCH_NOTES.md methodology):

weekly aggregation -> log transform -> STL seasonal adjustment ->

spike/outlier handling -> normalisation (logit share + per-capita).

We fit models on the SEASONALLY-ADJUSTED, UNSMOOTHED series; rolling means are

for the eye only (smoothing before change-point detection blurs break dates).

"""

from __future__ import annotations

import numpy as np

import pandas as pd

from statsmodels.tsa.seasonal import STL

def to_weekly(df_daily: pd.DataFrame) -> pd.DataFrame:

"""Sum to week-ending-Sunday. Removes the weekday CI cycle for free.

Drops the first/last partial weeks so every row is a full 7-day sum.

"""

wk = df_daily.resample("W-SUN").sum(min_count=7)

return wk.dropna(how="all")

def log_transform(df: pd.DataFrame) -> pd.DataFrame:

return np.log1p(df)

def deseasonalize(series: pd.Series, period: int = 52,

robust: bool = True) -> pd.Series:

"""Return STL seasonally-adjusted series (observed - seasonal) in the

same units as the input. Expects a regular weekly series."""

s = series.dropna()

if len(s) < 2 * period:

return s # too short to decompose

stl = STL(s, period=period, robust=robust).fit()

return s - stl.seasonal

def flag_spikes(resid: pd.Series, k: float = 3.5) -> pd.Series:

"""Robust MAD outlier flag on a residual series. True = anomalous

(e.g. release-week mirror stampede)."""

med = resid.median()

mad = (resid - med).abs().median()

if mad == 0:

return pd.Series(False, index=resid.index)

z = 0.6745 * (resid - med) / mad

return z.abs() > k

def winsorize_spikes(series: pd.Series, k: float = 3.5) -> pd.Series:

"""Cap robust-z outliers at the +/-k threshold (level series)."""

med = series.median()

mad = (series - med).abs().median()

if mad == 0:

return series

lo = med - k * mad / 0.6745

hi = med + k * mad / 0.6745

return series.clip(lo, hi)

def logit_share(df_levels: pd.DataFrame) -> pd.DataFrame:

"""Share-of-cohort in log-odds space. df_levels = wide weekly LEVELS

for the framework cohort (positive). Returns logit(share)."""

total = df_levels.sum(axis=1)

share = df_levels.div(total, axis=0).clip(1e-9, 1 - 1e-9)

return np.log(share / (1 - share))

def share_of_cohort(df_levels: pd.DataFrame) -> pd.DataFrame:

"""Plain fractional share (sums to 1 across columns)."""

total = df_levels.sum(axis=1)

return df_levels.div(total, axis=0)

def per_capita(df_levels: pd.DataFrame, denominator: pd.Series) -> pd.DataFrame:

"""Divide each framework by an ecosystem denominator (e.g. total npm

registry downloads) to detrend generic ecosystem growth."""

return df_levels.div(denominator, axis=0)

def index_to_baseline(df: pd.DataFrame, anchor: str) -> pd.DataFrame:

"""Rebase each series to 100 at the week nearest `anchor`."""

idx = df.index.get_indexer([pd.Timestamp(anchor)], method="nearest")[0]

base = df.iloc[idx]

return df.div(base) * 100.0

def dlog(df: pd.DataFrame) -> pd.DataFrame:

"""Week-over-week log growth (~ % growth)."""

return np.log(df).diff()

def rolling(df: pd.DataFrame, weeks: int = 4) -> pd.DataFrame:

"""Centred rolling mean -- VISUALISATION ONLY."""

return df.rolling(weeks, center=True, min_periods=1).mean()

def cagr(series: pd.Series, start: str, end: str) -> float:

"""Compound annual growth rate of a level series between two dates,

using the nearest available weeks."""

s = series.dropna()

i0 = s.index.get_indexer([pd.Timestamp(start)], method="nearest")[0]

i1 = s.index.get_indexer([pd.Timestamp(end)], method="nearest")[0]

v0, v1 = s.iloc[i0], s.iloc[i1]

years = (s.index[i1] - s.index[i0]).days / 365.25

if v0 <= 0 or years <= 0:

return float("nan")

return (v1 / v0) ** (1 / years) - 1.0