modeling_zhaocui/final/Q1/compare_models.py

from __future__ import annotations

import argparse
import csv
import hashlib
import importlib.metadata
import json
import math
import platform
import re
import subprocess
import sys
import time
import urllib.request
from dataclasses import dataclass
from pathlib import Path
from typing import Any, Iterable

import cv2
import matplotlib

matplotlib.use("Agg")
import matplotlib.pyplot as plt
import mediapipe as mp
import numpy as np
import torch
from openpyxl import load_workbook
from scipy.fft import dct, irfft, rfft
from scipy.spatial.transform import Rotation
from sklearn.linear_model import LogisticRegression, Ridge
from sklearn.metrics import accuracy_score, f1_score, mean_absolute_error
from sklearn.model_selection import GroupKFold
from transformers import AutoModel, AutoModelForCTC, AutoTokenizer


ROOT = Path(__file__).resolve().parents[2]
MATH_DIR = Path(__file__).resolve().parent
DATA_DIR = ROOT / "E棰樻暟鎹? / "闄勪欢1-鏁版嵁闆嗗師濮嬪妯℃€佹牱鏈? / "MOSEI鏁版嵁闆嗛儴鍒嗗師濮嬭棰?100鏉?
LABEL_FILE = DATA_DIR / "label-100.xlsx"
OUTPUT_DIR = MATH_DIR / "results" / "model_comparison"
CACHE_DIR = MATH_DIR / "cache" / "native"
FACE_MODEL_PATH = MATH_DIR / "cache" / "face_landmarker.task"

TEXT_MODEL_ID = "google-bert/bert-base-uncased"
SPEECH_MODEL_ID = "facebook/wav2vec2-base-960h"
SAMPLE_RATE = 16_000
FRAME_LENGTH = 400
FRAME_STEP = 160
N_FFT = 512
MEL_COUNT = 40
GRID_STEP_S = 0.1
VISION_RATE_HZ = 5.0
CTC_FRAME_STEP_S = 320 / SAMPLE_RATE
CTC_RECEPTIVE_SAMPLES = 400
SEED = 20260923
N_FOLDS = 5

ACTION_NAMES = (
    "browDownLeft", "browDownRight", "browInnerUp", "browOuterUpLeft",
    "browOuterUpRight", "eyeBlinkLeft", "eyeBlinkRight", "eyeSquintLeft",
    "eyeSquintRight", "eyeWideLeft", "eyeWideRight", "jawOpen",
    "mouthFrownLeft", "mouthFrownRight", "mouthPucker", "mouthSmileLeft",
    "mouthSmileRight",
)
VISION_NAMES = (
    *[f"mp_blendshape_{name}" for name in ACTION_NAMES],
    "pose_rotvec_x", "pose_rotvec_y", "pose_rotvec_z",
    "pose_translation_x", "pose_translation_y", "pose_translation_z",
    "gaze_left_x", "gaze_left_y", "gaze_left_z",
    "gaze_right_x", "gaze_right_y", "gaze_right_z",
    "eye_aperture_left", "eye_aperture_right", "mouth_aperture",
    "mouth_width_face_ratio", "brow_eye_distance_left", "brow_eye_distance_right",
)
AUDIO_NAMES = (
    *[f"logmel_{index:02d}" for index in range(40)],
    *[f"mfcc_{index:02d}" for index in range(13)],
    *[f"delta_mfcc_{index:02d}" for index in range(13)],
    "log_energy", "log_f0_hz", "voicing_strength", "spectral_centroid_hz",
    "spectral_bandwidth_hz", "spectral_flux", "zero_crossing_rate", "hnr_db",
)


@dataclass
class Sample:
    sample_id: str
    video_id: str
    clip_id: str
    duration_s: float
    sentiment: float
    polarity: int
    text_words: list[str]
    text_features: np.ndarray
    text_valid: np.ndarray
    hard_intervals: np.ndarray
    hard_valid: np.ndarray
    hard_quality: np.ndarray
    audio_times: np.ndarray
    audio_features: np.ndarray
    audio_observed: np.ndarray
    vision_times: np.ndarray
    vision_features: np.ndarray
    vision_observed: np.ndarray
    ctc_times: np.ndarray
    ctc_occupancy: np.ndarray
    speech_features: np.ndarray
    boundary_summary: list[dict[str, float]]
    video_sha256: str


@dataclass
class View:
    features: dict[str, np.ndarray]
    observed: dict[str, np.ndarray]
    coverage: dict[str, np.ndarray]
    edges: np.ndarray


def _sha256(path: Path) -> str:
    digest = hashlib.sha256()
    with path.open("rb") as stream:
        for chunk in iter(lambda: stream.read(1024 * 1024), b""):
            digest.update(chunk)
    return digest.hexdigest()


def _version(package: str) -> str:
    try:
        return importlib.metadata.version(package)
    except importlib.metadata.PackageNotFoundError:
        return "not-installed"


def _safe_name(value: str) -> str:
    return re.sub(r"[^A-Za-z0-9_.-]+", "_", value)


def _str_id(value: Any) -> str:
    if isinstance(value, float) and value.is_integer():
        return str(int(value))
    return str(value).strip()


def _read_labels(path: Path) -> list[dict[str, Any]]:
    workbook = load_workbook(path, read_only=True, data_only=True)
    sheet = workbook.active
    rows = sheet.iter_rows(values_only=True)
    header = [str(value).strip() if value is not None else "" for value in next(rows)]
    index = {name: position for position, name in enumerate(header)}
    required = {"video_id", "clip_id", "text", "label", "annotation"}
    missing = required - set(index)
    if missing:
        raise ValueError(f"label workbook lacks columns: {sorted(missing)}")
    records: list[dict[str, Any]] = []
    seen: set[str] = set()
    for row_number, row in enumerate(rows, start=2):
        video_id = _str_id(row[index["video_id"]])
        clip_id = _str_id(row[index["clip_id"]])
        sample_id = f"{video_id}/{clip_id}"
        if sample_id in seen:
            raise ValueError(f"duplicate source key in label workbook: {sample_id}")
        seen.add(sample_id)
        label = float(row[index["label"]])
        annotation = str(row[index["annotation"]]).strip().lower()
        polarity = 0 if label < 0 else 1 if label == 0 else 2
        expected = {"negative": 0, "neutral": 1, "positive": 2}.get(annotation)
        records.append({
            "row_number": row_number,
            "sample_id": sample_id,
            "video_id": video_id,
            "clip_id": clip_id,
            "text": str(row[index["text"]] or ""),
            "sentiment": label,
            "polarity": polarity,
            "annotation": annotation,
            "label_consistent": expected == polarity,
        })
    workbook.close()
    if len(records) != 100:
        raise ValueError(f"expected 100 label rows, received {len(records)}")
    return records


def _decode_audio(path: Path) -> np.ndarray:
    result = subprocess.run(
        [
            "ffmpeg", "-nostdin", "-v", "error", "-i", str(path), "-vn",
            "-ac", "1", "-ar", str(SAMPLE_RATE), "-f", "f32le", "pipe:1",
        ],
        check=True,
        capture_output=True,
    )
    audio = np.frombuffer(result.stdout, dtype="<f4").astype(np.float32, copy=True)
    if not len(audio):
        raise ValueError(f"no decoded audio samples in {path}")
    return np.nan_to_num(audio, nan=0.0, posinf=0.0, neginf=0.0)


def _duration(path: Path) -> float:
    result = subprocess.run(
        [
            "ffprobe", "-v", "error", "-show_entries", "format=duration",
            "-of", "default=noprint_wrappers=1:nokey=1", str(path),
        ],
        check=True,
        capture_output=True,
        text=True,
    )
    value = float(result.stdout.strip())
    if not math.isfinite(value) or value <= 0:
        raise ValueError(f"invalid duration for {path}: {value}")
    return value


def _mel_filterbank(sample_rate: int, n_fft: int, count: int) -> np.ndarray:
    def to_mel(freq: np.ndarray | float) -> np.ndarray | float:
        return 2595.0 * np.log10(1.0 + np.asarray(freq) / 700.0)

    def from_mel(mel: np.ndarray) -> np.ndarray:
        return 700.0 * (10.0 ** (mel / 2595.0) - 1.0)

    mel_points = np.linspace(to_mel(0.0), to_mel(sample_rate / 2), count + 2)
    hz_points = from_mel(mel_points)
    bins_hz = np.fft.rfftfreq(n_fft, 1 / sample_rate)
    bank = np.zeros((count, len(bins_hz)), dtype=np.float32)
    for band in range(count):
        left, center, right = hz_points[band:band + 3]
        bank[band] = np.maximum(
            0.0,
            np.minimum((bins_hz - left) / max(center - left, 1e-9),
                       (right - bins_hz) / max(right - center, 1e-9)),
        )
    return bank


def _audio_features(audio: np.ndarray, duration_s: float) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    starts = np.arange(0, max(len(audio), 1), FRAME_STEP, dtype=np.int64)
    frames = np.zeros((len(starts), FRAME_LENGTH), dtype=np.float32)
    for index, start in enumerate(starts):
        chunk = audio[start:start + FRAME_LENGTH]
        frames[index, :len(chunk)] = chunk
    window = np.hanning(FRAME_LENGTH).astype(np.float32)
    windowed = frames * window[None, :]
    spectrum = np.abs(rfft(windowed, n=N_FFT, axis=1)).astype(np.float32)
    power = spectrum ** 2
    mel = power @ _mel_filterbank(SAMPLE_RATE, N_FFT, MEL_COUNT).T
    log_mel = np.log(np.maximum(mel, 1e-8)).astype(np.float32)
    mfcc = dct(log_mel, type=2, norm="ortho", axis=1)[:, :13].astype(np.float32)
    delta = np.zeros_like(mfcc)
    if len(mfcc) > 1:
        for frame in range(len(mfcc)):
            left, right = max(0, frame - 2), min(len(mfcc), frame + 3)
            offsets = np.arange(left, right, dtype=np.float32) - frame
            denominator = float(np.square(offsets).sum())
            if denominator > 0:
                delta[frame] = (offsets[:, None] * mfcc[left:right]).sum(axis=0) / denominator

    frequencies = np.fft.rfftfreq(N_FFT, 1 / SAMPLE_RATE).astype(np.float32)
    magnitude_sum = np.maximum(spectrum.sum(axis=1), 1e-8)
    centroid = (spectrum * frequencies[None, :]).sum(axis=1) / magnitude_sum
    bandwidth = np.sqrt(
        (spectrum * (frequencies[None, :] - centroid[:, None]) ** 2).sum(axis=1)
        / magnitude_sum
    )
    normalized_spectrum = spectrum / magnitude_sum[:, None]
    flux = np.zeros(len(frames), dtype=np.float32)
    if len(frames) > 1:
        flux[1:] = np.square(normalized_spectrum[1:] - normalized_spectrum[:-1]).sum(axis=1)
    zcr = (np.diff(np.signbit(frames), axis=1) != 0).mean(axis=1).astype(np.float32)
    energy = np.mean(np.square(frames), axis=1)

    log_f0 = np.zeros(len(frames), dtype=np.float32)
    voicing = np.zeros(len(frames), dtype=np.float32)
    hnr = np.zeros(len(frames), dtype=np.float32)
    voiced = np.zeros(len(frames), dtype=np.bool_)
    min_lag, max_lag = int(SAMPLE_RATE / 400), int(SAMPLE_RATE / 60)
    for frame_index, frame in enumerate(windowed):
        if energy[frame_index] < 1e-7:
            continue
        autocorr = irfft(np.abs(rfft(frame, n=1024)) ** 2, n=1024)[:max_lag + 1]
        if autocorr[0] <= 1e-10:
            continue
        autocorr /= autocorr[0]
        region = autocorr[min_lag:max_lag + 1]
        if not len(region):
            continue
        lag = min_lag + int(np.argmax(region))
        strength = float(np.clip(autocorr[lag], 0.0, 1.0))
        voicing[frame_index] = strength
        voiced[frame_index] = strength >= 0.30
        if voiced[frame_index]:
            log_f0[frame_index] = math.log(SAMPLE_RATE / lag)
            hnr[frame_index] = 10.0 * math.log10(max(strength, 1e-5) / max(1.0 - strength, 1e-5))

    prosody = np.column_stack((
        np.log(np.maximum(energy, 1e-8)), log_f0, voicing, centroid, bandwidth,
        flux, zcr, hnr,
    )).astype(np.float32)
    features = np.column_stack((log_mel, mfcc, delta, prosody)).astype(np.float32)
    if features.shape[1] != 74:
        raise AssertionError(f"internal audio dimension error: {features.shape}")
    observed = np.isfinite(features)
    observed[:, 67] = voiced
    observed[:, 73] = voiced
    if len(frames) < 2:
        observed[:, 53:66] = False
    features = np.nan_to_num(features, nan=0.0, posinf=0.0, neginf=0.0)
    times = np.minimum((starts + FRAME_LENGTH / 2) / SAMPLE_RATE, duration_s).astype(np.float32)
    return times, features, observed


def _text_features(words: list[str], tokenizer: Any, model: Any, device: torch.device) -> tuple[np.ndarray, np.ndarray]:
    output = np.zeros((len(words), int(model.config.hidden_size)), dtype=np.float32)
    valid = np.zeros(len(words), dtype=np.bool_)
    if not words:
        return output, valid
    token_ids: list[int] = []
    word_ids: list[int] = []
    for word_index, word in enumerate(words):
        word_tokens = tokenizer(word, add_special_tokens=False)["input_ids"]
        token_ids.extend(int(value) for value in word_tokens)
        word_ids.extend([word_index] * len(word_tokens))
    token_sums = np.zeros((len(token_ids), output.shape[1]), dtype=np.float64)
    token_weights = np.zeros(len(token_ids), dtype=np.float64)
    window_size, stride = 510, 384
    offsets = range(0, max(len(token_ids), 1), stride)
    with torch.inference_mode():
        for offset in offsets:
            end = min(len(token_ids), offset + window_size)
            if offset >= end:
                break
            content = token_ids[offset:end]
            special_content = [int(tokenizer.cls_token_id), *content, int(tokenizer.sep_token_id)]
            encoded = {
                "input_ids": torch.tensor([special_content], dtype=torch.long, device=device),
                "attention_mask": torch.ones((1, len(special_content)), dtype=torch.long, device=device),
            }
            hidden_states = model(**encoded, output_hidden_states=True).hidden_states
            hidden = torch.stack(tuple(hidden_states[-4:]), dim=0).mean(dim=0)[0, 1:1 + len(content)]
            values = hidden.float().cpu().numpy()
            for local_index, value in enumerate(values):
                global_index = offset + local_index
                edge_weight = float(min(local_index + 1, len(content) - local_index))
                token_sums[global_index] += edge_weight * value
                token_weights[global_index] += edge_weight
    accum = np.zeros_like(output, dtype=np.float64)
    counts = np.zeros(len(words), dtype=np.float64)
    for token_index, word_index in enumerate(word_ids):
        if word_index is None or not (0 <= word_index < len(words)) or token_weights[token_index] <= 0:
            continue
        accum[word_index] += token_sums[token_index] / token_weights[token_index]
        counts[word_index] += 1.0
    valid = counts > 0
    output[valid] = (accum[valid] / counts[valid, None]).astype(np.float32)
    return output, valid


def _ctc_targets(words: list[str], tokenizer: Any) -> tuple[list[int], list[list[int]]]:
    vocab = tokenizer.get_vocab()
    delimiter_id = int(tokenizer.convert_tokens_to_ids(tokenizer.word_delimiter_token or "|"))
    unknown_id = int(tokenizer.unk_token_id)
    targets: list[int] = []
    per_word: list[list[int]] = [[] for _ in words]
    for word_index, raw_word in enumerate(words):
        if word_index:
            targets.append(delimiter_id)
        normalized = re.sub(r"[^a-z']", "", raw_word.lower().replace("鈥?, "'"))
        for char in normalized:
            per_word[word_index].append(len(targets))
            targets.append(int(vocab.get(char, unknown_id)))
    return targets, per_word


def _ctc_path_and_posterior(log_probs: np.ndarray, targets: list[int], blank_id: int) -> tuple[np.ndarray | None, np.ndarray | None, float]:
    if not targets or log_probs.ndim != 2:
        return None, None, -math.inf
    target = np.asarray(targets, dtype=np.int64)
    states = np.full(2 * len(target) + 1, blank_id, dtype=np.int64)
    states[1::2] = target
    frames, state_count = log_probs.shape[0], len(states)
    if frames == 0 or frames < len(target):
        return None, None, -math.inf
    alpha = np.full((frames, state_count), -np.inf, dtype=np.float64)
    alpha[0, 0] = float(log_probs[0, blank_id])
    alpha[0, 1] = float(log_probs[0, states[1]])
    skip_ok = np.zeros(state_count, dtype=np.bool_)
    if state_count > 2:
        indices = np.arange(state_count)
        skip_ok[2:] = (states[2:] != blank_id) & (states[2:] != states[:-2])
    for frame in range(1, frames):
        prev = alpha[frame - 1]
        incoming = prev.copy()
        incoming[1:] = np.logaddexp(incoming[1:], prev[:-1])
        valid_skip = np.flatnonzero(skip_ok)
        if len(valid_skip):
            incoming[valid_skip] = np.logaddexp(incoming[valid_skip], prev[valid_skip - 2])
        alpha[frame] = incoming + log_probs[frame, states]
    log_z = float(np.logaddexp(alpha[-1, -1], alpha[-1, -2]))
    if not math.isfinite(log_z):
        return None, None, -math.inf

    beta = np.full((frames, state_count), -np.inf, dtype=np.float64)
    beta[-1, -1] = 0.0
    beta[-1, -2] = 0.0
    for frame in range(frames - 2, -1, -1):
        nxt = beta[frame + 1] + log_probs[frame + 1, states]
        # A CTC state may stay active for any number of frames. Start with the
        # same-state transition, then add transitions to the next states.
        values = nxt.copy()
        values[:-1] = np.logaddexp(values[:-1], nxt[1:])
        if state_count > 2:
            dest = np.flatnonzero(skip_ok)
            values[dest - 2] = np.logaddexp(values[dest - 2], nxt[dest])
        beta[frame] = values
    posterior = np.exp(np.clip(alpha + beta - log_z, -745.0, 0.0)).astype(np.float32)

    back = np.zeros((frames, state_count), dtype=np.uint8)
    previous = np.full(state_count, -np.inf, dtype=np.float64)
    previous[0] = float(log_probs[0, blank_id])
    previous[1] = float(log_probs[0, states[1]])
    for frame in range(1, frames):
        stay = previous
        one = np.full(state_count, -np.inf, dtype=np.float64)
        one[1:] = previous[:-1]
        two = np.full(state_count, -np.inf, dtype=np.float64)
        two[skip_ok] = previous[np.flatnonzero(skip_ok) - 2]
        candidates = np.stack((stay, one, two), axis=0)
        choice = candidates.argmax(axis=0).astype(np.uint8)
        previous = candidates[choice, np.arange(state_count)] + log_probs[frame, states]
        back[frame] = choice
    state = state_count - 1 if previous[-1] >= previous[-2] else state_count - 2
    path = np.empty(frames, dtype=np.int32)
    path[-1] = state
    for frame in range(frames - 1, 0, -1):
        state -= int(back[frame, state])
        path[frame - 1] = state
    return path, posterior, log_z


def _entry_distribution(log_probs: np.ndarray, states: np.ndarray, alpha: np.ndarray, beta: np.ndarray, log_z: float, state: int) -> np.ndarray:
    frames = len(log_probs)
    mass = np.zeros(frames, dtype=np.float64)
    if state == 1:
        mass[0] = math.exp(float(log_probs[0, states[state]] + beta[0, state] - log_z))
    predecessors = [state - 1]
    if state >= 2 and states[state] != states[state - 2]:
        predecessors.append(state - 2)
    for frame in range(1, frames):
        terms = [alpha[frame - 1, prev] for prev in predecessors if prev >= 0]
        if terms:
            log_mass = float(np.logaddexp.reduce(terms) + log_probs[frame, states[state]] + beta[frame, state] - log_z)
            mass[frame] = math.exp(max(-745.0, min(log_mass, 0.0)))
    return mass


def _exit_distribution(log_probs: np.ndarray, states: np.ndarray, alpha: np.ndarray, beta: np.ndarray, log_z: float, state: int) -> np.ndarray:
    frames = len(log_probs)
    mass = np.zeros(frames + 1, dtype=np.float64)
    successors = [state + 1]
    if state + 2 < len(states) and states[state] != states[state + 2]:
        successors.append(state + 2)
    for frame in range(1, frames):
        terms = [
            alpha[frame - 1, state] + log_probs[frame, states[nxt]] + beta[frame, nxt]
            for nxt in successors if nxt < len(states)
        ]
        if terms:
            mass[frame] = math.exp(max(-745.0, min(float(np.logaddexp.reduce(terms) - log_z), 0.0)))
    if state == len(states) - 2:
        mass[frames] = math.exp(max(-745.0, min(float(alpha[-1, state] - log_z), 0.0)))
    return mass


def _weighted_quantile(values: np.ndarray, weights: np.ndarray, quantile: float) -> float:
    if not len(values) or float(weights.sum()) <= 0:
        return math.nan
    order = np.argsort(values)
    values, weights = values[order], weights[order]
    cumulative = np.cumsum(weights) / weights.sum()
    return float(values[min(int(np.searchsorted(cumulative, quantile, side="left")), len(values) - 1)])


def _ctc_boundaries(
    log_probs: np.ndarray, targets: list[int], word_targets: list[list[int]],
    blank_id: int, duration_s: float, center_s: float,
) -> tuple[np.ndarray, list[dict[str, float]], np.ndarray]:
    path, posterior, log_z = _ctc_path_and_posterior(log_probs, targets, blank_id)
    occupancy = np.zeros((len(log_probs), len(word_targets)), dtype=np.float32)
    summaries: list[dict[str, float]] = [
        {"start_mean_s": math.nan, "end_mean_s": math.nan,
         "start_p05_s": math.nan, "start_p95_s": math.nan,
         "end_p05_s": math.nan, "end_p95_s": math.nan,
         "start_width90_s": math.nan, "end_width90_s": math.nan}
        for _ in word_targets
    ]
    if path is None or posterior is None:
        return np.empty((0, len(word_targets)), dtype=np.float32), summaries, path if path is not None else np.empty((0,), dtype=np.int32)
    state_ids = np.full(2 * len(targets) + 1, blank_id, dtype=np.int64)
    state_ids[1::2] = np.asarray(targets, dtype=np.int64)
    # Recompute the dynamic-programming tables for boundary transition marginals.
    frames, count = len(log_probs), len(state_ids)
    alpha = np.full((frames, count), -np.inf, dtype=np.float64)
    alpha[0, 0] = log_probs[0, blank_id]
    alpha[0, 1] = log_probs[0, state_ids[1]]
    skip = np.zeros(count, dtype=np.bool_)
    if count > 2:
        skip[2:] = (state_ids[2:] != blank_id) & (state_ids[2:] != state_ids[:-2])
    for frame in range(1, frames):
        incoming = alpha[frame - 1].copy()
        incoming[1:] = np.logaddexp(incoming[1:], alpha[frame - 1, :-1])
        skip_states = np.flatnonzero(skip)
        if len(skip_states):
            incoming[skip_states] = np.logaddexp(incoming[skip_states], alpha[frame - 1, skip_states - 2])
        alpha[frame] = incoming + log_probs[frame, state_ids]
    beta = np.full((frames, count), -np.inf, dtype=np.float64)
    beta[-1, -1] = beta[-1, -2] = 0.0
    for frame in range(frames - 2, -1, -1):
        nxt = beta[frame + 1] + log_probs[frame + 1, state_ids]
        current = nxt.copy()
        current[:-1] = np.logaddexp(current[:-1], nxt[1:])
        skip_states = np.flatnonzero(skip)
        if len(skip_states):
            current[skip_states - 2] = np.logaddexp(current[skip_states - 2], nxt[skip_states])
        beta[frame] = current
    frame_step = CTC_FRAME_STEP_S
    first_edge = max(0.0, center_s - frame_step / 2)
    frame_edges = np.clip(first_edge + np.arange(frames + 1, dtype=np.float64) * frame_step, 0.0, duration_s)
    for word_index, token_indices in enumerate(word_targets):
        if not token_indices:
            continue
        char_states = np.asarray([2 * target_index + 1 for target_index in token_indices], dtype=np.int64)
        occupancy[:, word_index] = posterior[:, char_states].sum(axis=1)
        start_mass = _entry_distribution(log_probs, state_ids, alpha, beta, log_z, int(char_states[0]))
        end_mass = _exit_distribution(log_probs, state_ids, alpha, beta, log_z, int(char_states[-1]))
        start_total = float(start_mass.sum())
        end_total = float(end_mass.sum())
        if not math.isclose(start_total, 1.0, rel_tol=2e-4, abs_tol=2e-4):
            raise FloatingPointError(f"CTC start-boundary marginal sums to {start_total:.6g}, expected 1")
        if not math.isclose(end_total, 1.0, rel_tol=2e-4, abs_tol=2e-4):
            raise FloatingPointError(f"CTC end-boundary marginal sums to {end_total:.6g}, expected 1")
        start_times = frame_edges[:-1]
        end_times = frame_edges
        start_p05 = _weighted_quantile(start_times, start_mass, 0.05)
        start_p95 = _weighted_quantile(start_times, start_mass, 0.95)
        end_p05 = _weighted_quantile(end_times, end_mass, 0.05)
        end_p95 = _weighted_quantile(end_times, end_mass, 0.95)
        summaries[word_index] = {
            "start_mean_s": float(np.sum(start_mass * start_times) / start_total),
            "end_mean_s": float(np.sum(end_mass * end_times) / end_total),
            "start_p05_s": start_p05,
            "start_p95_s": start_p95,
            "end_p05_s": end_p05,
            "end_p95_s": end_p95,
            "start_width90_s": start_p95 - start_p05,
            "end_width90_s": end_p95 - end_p05,
        }
    return occupancy, summaries, path


def _posterior_hard_intervals(
    log_probs: np.ndarray,
    path: np.ndarray | None,
    targets: list[int],
    word_targets: list[list[int]],
    duration_s: float,
    frame_step: float,
    center_s: float,
) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    intervals = np.zeros((len(word_targets), 2), dtype=np.float32)
    valid = np.zeros(len(word_targets), dtype=np.bool_)
    quality = np.zeros(len(word_targets), dtype=np.float32)
    if path is None or not len(path):
        return intervals, valid, quality
    for word_index, target_indices in enumerate(word_targets):
        if not target_indices:
            continue
        char_states = np.asarray([2 * index + 1 for index in target_indices], dtype=np.int32)
        frame_indices = np.flatnonzero(np.isin(path, char_states))
        if not len(frame_indices):
            continue
        first, last = int(frame_indices[0]), int(frame_indices[-1])
        start = max(0.0, first * frame_step + center_s - frame_step / 2)
        end = min(duration_s, (last + 1) * frame_step + center_s - frame_step / 2)
        if end <= start:
            continue
        intervals[word_index] = (start, end)
        valid[word_index] = True
        scores = []
        for target_index in target_indices:
            selected = np.flatnonzero(path == 2 * target_index + 1)
            if len(selected):
                scores.extend(log_probs[selected, targets[target_index]].tolist())
        quality[word_index] = float(np.exp(np.mean(scores))) if scores else 0.0
    return intervals, valid, quality


def _landmark_points(landmarks: Iterable[Any]) -> np.ndarray:
    return np.asarray([[item.x, item.y, item.z] for item in landmarks], dtype=np.float32)


def _distance(points: np.ndarray, left: int, right: int) -> float:
    return float(np.linalg.norm(points[left] - points[right]))


def _gaze_proxy(points: np.ndarray, left_eye: bool) -> np.ndarray | None:
    if points.shape[0] < 478:
        return None
    if left_eye:
        corner_a, corner_b, upper, lower = 33, 133, 159, 145
        iris_ids = np.arange(468, 473)
    else:
        corner_a, corner_b, upper, lower = 362, 263, 386, 374
        iris_ids = np.arange(473, 478)
    horizontal = points[corner_b] - points[corner_a]
    width = float(np.linalg.norm(horizontal))
    vertical = points[upper] - points[lower]
    height = float(np.linalg.norm(vertical))
    if width <= 1e-6 or height <= 1e-6:
        return None
    u = horizontal / width
    v = vertical - float(np.dot(vertical, u)) * u
    v_norm = float(np.linalg.norm(v))
    if v_norm <= 1e-6:
        return None
    v /= v_norm
    eye_center = (points[corner_a] + points[corner_b]) / 2
    iris_center = points[iris_ids].mean(axis=0)
    delta = iris_center - eye_center
    gx = float(np.dot(delta, u) / width)
    gy = float(np.dot(delta, v) / width)
    # Face-relative iris offset plus a fixed forward component; this is a gaze proxy.
    direction = np.asarray([-gx, -gy, -1.0], dtype=np.float32)
    norm = float(np.linalg.norm(direction))
    if norm <= 1e-8:
        return None
    return direction / norm


def _pose_from_matrix(matrix_obj: Any) -> tuple[np.ndarray, np.ndarray] | None:
    raw = np.asarray(matrix_obj.data, dtype=np.float64)
    rows = int(getattr(matrix_obj, "rows", 4) or 4)
    cols = int(getattr(matrix_obj, "cols", 4) or 4)
    if raw.size != rows * cols or rows < 3 or cols < 4:
        return None
    matrix = raw.reshape(rows, cols)
    rotation_raw = matrix[:3, :3]
    if not np.isfinite(rotation_raw).all() or not np.isfinite(matrix[:3, 3]).all():
        return None
    u, _, vh = np.linalg.svd(rotation_raw)
    rotation_matrix = u @ vh
    if np.linalg.det(rotation_matrix) < 0:
        u[:, -1] *= -1
        rotation_matrix = u @ vh
    return Rotation.from_matrix(rotation_matrix).as_rotvec().astype(np.float32), matrix[:3, 3].astype(np.float32)


def _vision_frame_features(landmarks: list[Any], blendshapes: list[Any], matrix_obj: Any | None) -> tuple[np.ndarray, np.ndarray]:
    features = np.zeros(35, dtype=np.float32)
    observed = np.zeros(35, dtype=np.bool_)
    shape_lookup = {category.category_name: float(category.score) for category in blendshapes}
    for index, name in enumerate(ACTION_NAMES):
        if name in shape_lookup and math.isfinite(shape_lookup[name]):
            features[index] = shape_lookup[name]
            observed[index] = True

    if matrix_obj is not None:
        pose = _pose_from_matrix(matrix_obj)
        if pose is not None:
            features[17:20] = pose[0]
            features[20:23] = pose[1]
            observed[17:23] = True

    points = _landmark_points(landmarks)
    for eye_index, is_left in enumerate((True, False)):
        gaze = _gaze_proxy(points, is_left)
        if gaze is not None:
            start = 23 + eye_index * 3
            features[start:start + 3] = gaze
            observed[start:start + 3] = True
    if points.shape[0] >= 455:
        face_width = _distance(points, 234, 454)
        left_eye_width = _distance(points, 33, 133)
        right_eye_width = _distance(points, 362, 263)
        mouth_width = _distance(points, 61, 291)
        if min(face_width, left_eye_width, right_eye_width, mouth_width) > 1e-6:
            ratios = np.asarray((
                _distance(points, 159, 145) / left_eye_width,
                _distance(points, 386, 374) / right_eye_width,
                _distance(points, 13, 14) / mouth_width,
                mouth_width / face_width,
                _distance(points, 105, 159) / left_eye_width,
                _distance(points, 334, 386) / right_eye_width,
            ), dtype=np.float32)
            features[29:35] = ratios
            observed[29:35] = np.isfinite(ratios)
    return features, observed


def _ensure_face_model() -> str:
    FACE_MODEL_PATH.parent.mkdir(parents=True, exist_ok=True)
    if not FACE_MODEL_PATH.is_file() or FACE_MODEL_PATH.stat().st_size < 1_000_000:
        url = (
            "https://storage.googleapis.com/mediapipe-models/face_landmarker/"
            "face_landmarker/float16/latest/face_landmarker.task"
        )
        temporary = FACE_MODEL_PATH.with_suffix(".task.tmp")
        urllib.request.urlretrieve(url, temporary)
        if temporary.stat().st_size < 1_000_000:
            temporary.unlink(missing_ok=True)
            raise RuntimeError("downloaded Face Landmarker model is unexpectedly small")
        temporary.replace(FACE_MODEL_PATH)
    return _sha256(FACE_MODEL_PATH)


def _vision_features(path: Path, duration_s: float, face_model_path: Path) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    options = mp.tasks.vision.FaceLandmarkerOptions(
        base_options=mp.tasks.BaseOptions(model_asset_path=str(face_model_path)),
        running_mode=mp.tasks.vision.RunningMode.VIDEO,
        num_faces=1,
        output_face_blendshapes=True,
        output_facial_transformation_matrixes=True,
    )
    capture = cv2.VideoCapture(str(path))
    if not capture.isOpened():
        raise RuntimeError(f"OpenCV cannot decode video: {path}")
    fps = float(capture.get(cv2.CAP_PROP_FPS) or 0.0)
    if not math.isfinite(fps) or fps <= 0:
        fps = 30.0
    timestamps: list[float] = []
    values: list[np.ndarray] = []
    masks: list[np.ndarray] = []
    next_time = 0.0
    frame_index = 0
    last_timestamp_ms = -1
    with mp.tasks.vision.FaceLandmarker.create_from_options(options) as landmarker:
        while True:
            ok, bgr = capture.read()
            if not ok:
                break
            reported_ms = float(capture.get(cv2.CAP_PROP_POS_MSEC))
            timestamp = reported_ms / 1000.0
            previous_timestamp = timestamps[-1] if timestamps else -1.0
            if not math.isfinite(timestamp) or (frame_index > 0 and timestamp <= previous_timestamp):
                timestamp = frame_index / fps
            frame_index += 1
            if timestamp + 1e-6 < next_time:
                continue
            while next_time <= timestamp:
                next_time += 1.0 / VISION_RATE_HZ
            timestamp = min(max(timestamp, 0.0), duration_s)
            timestamp_ms = max(last_timestamp_ms + 1, int(round(timestamp * 1000)))
            last_timestamp_ms = timestamp_ms
            rgb = cv2.cvtColor(bgr, cv2.COLOR_BGR2RGB)
            image = mp.Image(image_format=mp.ImageFormat.SRGB, data=np.ascontiguousarray(rgb))
            result = landmarker.detect_for_video(image, timestamp_ms)
            timestamps.append(timestamp)
            faces = result.face_landmarks
            if not faces:
                values.append(np.zeros(35, dtype=np.float32))
                masks.append(np.zeros(35, dtype=np.bool_))
                continue
            blendshapes = result.face_blendshapes[0] if result.face_blendshapes else []
            matrix = result.facial_transformation_matrixes[0] if result.facial_transformation_matrixes else None
            vector, mask = _vision_frame_features(faces[0], blendshapes, matrix)
            values.append(vector)
            masks.append(mask)
    capture.release()
    if not timestamps:
        return np.empty(0, np.float32), np.empty((0, 35), np.float32), np.empty((0, 35), np.bool_)
    return np.asarray(timestamps, np.float32), np.stack(values), np.stack(masks)


def _load_models(device: torch.device) -> tuple[Any, Any, Any, Any, dict[str, Any]]:
    text_tokenizer = AutoTokenizer.from_pretrained(TEXT_MODEL_ID, use_fast=True)
    text_model = AutoModel.from_pretrained(TEXT_MODEL_ID, output_hidden_states=True).to(device).eval()
    speech_tokenizer = AutoTokenizer.from_pretrained(SPEECH_MODEL_ID)
    speech_model = AutoModelForCTC.from_pretrained(SPEECH_MODEL_ID).to(device).eval()
    config = speech_model.config
    stride_samples = int(np.prod(config.conv_stride))
    receptive = 1
    jump = 1
    for kernel, stride in zip(config.conv_kernel, config.conv_stride):
        receptive += (int(kernel) - 1) * jump
        jump *= int(stride)
    info = {
        "text_id": TEXT_MODEL_ID,
        "text_revision": getattr(text_model.config, "_commit_hash", None),
        "speech_id": SPEECH_MODEL_ID,
        "speech_revision": getattr(config, "_commit_hash", None),
        "speech_stride_samples": stride_samples,
        "speech_receptive_field_samples": receptive,
    }
    return text_tokenizer, text_model, speech_tokenizer, speech_model, info


def _extract_native_sample(
    record: dict[str, Any],
    models: tuple[Any, Any, Any, Any, dict[str, Any]],
    device: torch.device,
    face_model_path: Path,
    video_sha256: str,
) -> Sample:
    text_tokenizer, text_model, ctc_tokenizer, speech_model, model_info = models
    video_path = DATA_DIR / record["video_id"] / f"{record['clip_id']}.mp4"
    duration_s = _duration(video_path)
    words = record["text"].split()
    text_vectors, text_valid = _text_features(words, text_tokenizer, text_model, device)
    waveform = _decode_audio(video_path)
    audio_times, audio_vectors, audio_observed = _audio_features(waveform, duration_s)
    vision_times, vision_vectors, vision_observed = _vision_features(video_path, duration_s, face_model_path)
    targets, word_targets = _ctc_targets(words, ctc_tokenizer)
    blank_id = int(ctc_tokenizer.pad_token_id)
    if targets:
        with torch.inference_mode():
            input_values = torch.from_numpy(waveform).to(device).unsqueeze(0)
            output = speech_model(input_values=input_values, output_hidden_states=True)
            logits = output.logits[0].float()
            log_probs = torch.log_softmax(logits, dim=-1).cpu().numpy().astype(np.float32)
            deep_states = torch.stack(tuple(output.hidden_states[-4:]), dim=0).mean(dim=0)[0].float().cpu().numpy()
    else:
        log_probs = np.empty((0, len(ctc_tokenizer)), dtype=np.float32)
        deep_states = np.empty((0, int(speech_model.config.hidden_size)), dtype=np.float32)
    stride_samples = int(model_info["speech_stride_samples"])
    frame_step = stride_samples / SAMPLE_RATE
    center_s = float(model_info["speech_receptive_field_samples"]) / (2 * SAMPLE_RATE)
    if len(log_probs) and targets:
        occupancy, boundary_summary, path = _ctc_boundaries(
            log_probs, targets, word_targets, blank_id, duration_s, center_s
        )
        hard_intervals, hard_valid, hard_quality = _posterior_hard_intervals(
            log_probs, path, targets, word_targets, duration_s, frame_step, center_s
        )
    else:
        occupancy = np.zeros((len(log_probs), len(words)), dtype=np.float32)
        boundary_summary = [
        {"start_mean_s": math.nan, "end_mean_s": math.nan,
         "start_p05_s": math.nan, "start_p95_s": math.nan,
         "end_p05_s": math.nan, "end_p95_s": math.nan,
         "start_width90_s": math.nan, "end_width90_s": math.nan}
            for _ in words
        ]
        hard_intervals = np.zeros((len(words), 2), dtype=np.float32)
        hard_valid = np.zeros(len(words), dtype=np.bool_)
        hard_quality = np.zeros(len(words), dtype=np.float32)
    ctc_frame_times = (np.arange(len(log_probs), dtype=np.float32) * frame_step + center_s).clip(0.0, duration_s)
    return Sample(
        sample_id=record["sample_id"],
        video_id=record["video_id"],
        clip_id=record["clip_id"],
        duration_s=duration_s,
        sentiment=record["sentiment"],
        polarity=record["polarity"],
        text_words=words,
        text_features=text_vectors.astype(np.float16),
        text_valid=text_valid,
        hard_intervals=hard_intervals,
        hard_valid=hard_valid,
        hard_quality=hard_quality,
        audio_times=audio_times,
        audio_features=audio_vectors.astype(np.float16),
        audio_observed=audio_observed,
        vision_times=vision_times,
        vision_features=vision_vectors.astype(np.float16),
        vision_observed=vision_observed,
        ctc_times=ctc_frame_times,
        ctc_occupancy=occupancy.astype(np.float16),
        speech_features=deep_states.astype(np.float16),
        boundary_summary=boundary_summary,
        video_sha256=video_sha256,
    )


def _save_cache(path: Path, sample: Sample, cache_schema: str) -> None:
    path.parent.mkdir(parents=True, exist_ok=True)
    temporary = path.with_suffix(".npz.tmp")
    with temporary.open("wb") as stream:
        np.savez_compressed(
            stream,
            cache_schema=np.asarray(cache_schema),
            sample_id=np.asarray(sample.sample_id),
            video_id=np.asarray(sample.video_id),
            clip_id=np.asarray(sample.clip_id),
            duration_s=np.asarray(sample.duration_s, dtype=np.float32),
            sentiment=np.asarray(sample.sentiment, dtype=np.float32),
            polarity=np.asarray(sample.polarity, dtype=np.int8),
            text_words=np.asarray(sample.text_words, dtype="U160"),
            text_features=sample.text_features.astype(np.float16),
            text_valid=sample.text_valid,
            hard_intervals=sample.hard_intervals.astype(np.float32),
            hard_valid=sample.hard_valid,
            hard_quality=sample.hard_quality.astype(np.float32),
            audio_times=sample.audio_times.astype(np.float32),
            audio_features=sample.audio_features.astype(np.float16),
            audio_observed=sample.audio_observed,
            vision_times=sample.vision_times.astype(np.float32),
            vision_features=sample.vision_features.astype(np.float16),
            vision_observed=sample.vision_observed,
            ctc_times=sample.ctc_times.astype(np.float32),
            ctc_occupancy=sample.ctc_occupancy.astype(np.float16),
            speech_features=sample.speech_features.astype(np.float16),
            boundary_summary=json.dumps(sample.boundary_summary, ensure_ascii=False),
            video_sha256=np.asarray(sample.video_sha256),
        )
    temporary.replace(path)


def _load_cache(path: Path, record: dict[str, Any], source_hash: str, cache_schema: str) -> Sample:
    with np.load(path, allow_pickle=False) as data:
        if str(data["cache_schema"].item()) != cache_schema:
            raise ValueError("feature cache schema mismatch")
        if str(data["sample_id"].item()) != record["sample_id"]:
            raise ValueError("feature cache sample key mismatch")
        if str(data["video_sha256"].item()) != source_hash:
            raise ValueError("feature cache source hash mismatch")
        return Sample(
            sample_id=record["sample_id"], video_id=record["video_id"], clip_id=record["clip_id"],
            duration_s=float(data["duration_s"]), sentiment=float(data["sentiment"]), polarity=int(data["polarity"]),
            text_words=data["text_words"].astype(str).tolist(), text_features=np.asarray(data["text_features"], dtype=np.float16),
            text_valid=np.asarray(data["text_valid"], dtype=np.bool_), hard_intervals=np.asarray(data["hard_intervals"], dtype=np.float32),
            hard_valid=np.asarray(data["hard_valid"], dtype=np.bool_), hard_quality=np.asarray(data["hard_quality"], dtype=np.float32),
            audio_times=np.asarray(data["audio_times"], dtype=np.float32), audio_features=np.asarray(data["audio_features"], dtype=np.float16),
            audio_observed=np.asarray(data["audio_observed"], dtype=np.bool_), vision_times=np.asarray(data["vision_times"], dtype=np.float32),
            vision_features=np.asarray(data["vision_features"], dtype=np.float16), vision_observed=np.asarray(data["vision_observed"], dtype=np.bool_),
            ctc_times=np.asarray(data["ctc_times"], dtype=np.float32), ctc_occupancy=np.asarray(data["ctc_occupancy"], dtype=np.float16),
            speech_features=np.asarray(data["speech_features"], dtype=np.float16),
            boundary_summary=json.loads(str(data["boundary_summary"].item())), video_sha256=source_hash,
        )


def _voronoi_intervals(times: np.ndarray, duration_s: float) -> np.ndarray:
    if not len(times):
        return np.empty((0, 2), dtype=np.float32)
    centers = np.clip(np.asarray(times, dtype=np.float64), 0.0, duration_s)
    if np.any(np.diff(centers) < -1e-7):
        order = np.argsort(centers)
        centers = centers[order]
    else:
        order = np.arange(len(centers))
    # Duplicate PTS values are coalesced to a tiny increasing interval for assignment.
    for index in range(1, len(centers)):
        if centers[index] <= centers[index - 1]:
            centers[index] = min(duration_s, centers[index - 1] + 1e-6)
    edges = np.empty(len(centers) + 1, dtype=np.float64)
    edges[0], edges[-1] = 0.0, duration_s
    if len(centers) > 1:
        edges[1:-1] = (centers[:-1] + centers[1:]) / 2
    edges = np.maximum.accumulate(np.clip(edges, 0.0, duration_s))
    intervals = np.column_stack((edges[:-1], edges[1:])).astype(np.float32)
    if len(order) != len(intervals) or not np.array_equal(order, np.arange(len(order))):
        restored = np.zeros_like(intervals)
        restored[order] = intervals
        return restored
    return intervals


def _grid_edges(duration_s: float) -> np.ndarray:
    count = max(1, int(math.ceil(duration_s / GRID_STEP_S)))
    edges = np.arange(count + 1, dtype=np.float64) * GRID_STEP_S
    edges[-1] = duration_s
    return edges


def _project_rows(
    values: np.ndarray,
    observed: np.ndarray,
    source_intervals: np.ndarray,
    qualities: np.ndarray,
    target_edges: np.ndarray,
) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    dim = values.shape[1] if values.ndim == 2 else 0
    count = len(target_edges) - 1
    output = np.zeros((count, dim), dtype=np.float32)
    mask = np.zeros((count, dim), dtype=np.bool_)
    coverage = np.zeros((count, dim), dtype=np.float32)
    if not dim or not len(source_intervals):
        return output, mask, coverage
    intervals = np.asarray(source_intervals, dtype=np.float64)
    qualities = np.maximum(np.asarray(qualities, dtype=np.float64), 0.0)
    values64 = np.asarray(values, dtype=np.float64)
    for target_index, (left, right) in enumerate(zip(target_edges[:-1], target_edges[1:])):
        if right <= left:
            continue
        overlap = np.maximum(
            0.0,
            np.minimum(intervals[:, 1], right) - np.maximum(intervals[:, 0], left),
        ) * qualities
        candidate = np.flatnonzero(overlap > 0)
        if not len(candidate):
            continue
        local_mask = observed[candidate]
        local_weight = overlap[candidate, None] * local_mask
        denominator = local_weight.sum(axis=0)
        physical_coverage = (overlap[candidate, None] * local_mask).sum(axis=0)
        good = denominator > 0
        if good.any():
            output[target_index, good] = (
                (values64[candidate] * local_weight).sum(axis=0)[good] / denominator[good]
            ).astype(np.float32)
            mask[target_index, good] = True
            coverage[target_index, good] = np.minimum(1.0, physical_coverage[good] / (right - left)).astype(np.float32)
    return output, mask, coverage


def _so3_weighted_mean(rotvecs: np.ndarray, weights: np.ndarray) -> np.ndarray | None:
    positive = weights > 0
    if not positive.any():
        return None
    rotations = Rotation.from_rotvec(np.asarray(rotvecs[positive], dtype=np.float64))
    local_weights = np.asarray(weights[positive], dtype=np.float64)
    local_weights /= local_weights.sum()
    matrices = rotations.as_matrix()
    mean_matrix = np.einsum("n,nij->ij", local_weights, matrices)
    u, _, vh = np.linalg.svd(mean_matrix)
    mean_matrix = u @ vh
    if np.linalg.det(mean_matrix) < 0:
        u[:, -1] *= -1
        mean_matrix = u @ vh
    mean_rotation = Rotation.from_matrix(mean_matrix)
    for _ in range(30):
        residual = (mean_rotation.inv() * rotations).as_rotvec()
        delta = np.einsum("n,nd->d", local_weights, residual)
        if np.linalg.norm(delta) < 1e-8:
            break
        mean_rotation = mean_rotation * Rotation.from_rotvec(delta)
    return mean_rotation.as_rotvec().astype(np.float32)


def _project_vision_geometry(sample: Sample, target_edges: np.ndarray, geometry_aware: bool) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    intervals = _voronoi_intervals(sample.vision_times, sample.duration_s)
    quality = np.ones(len(intervals), dtype=np.float32)
    values, masks, coverage = _project_rows(
        sample.vision_features.astype(np.float32), sample.vision_observed, intervals, quality, target_edges
    )
    if not geometry_aware or not len(intervals):
        return values, masks, coverage
    raw = sample.vision_features.astype(np.float32)
    for target_index, (left, right) in enumerate(zip(target_edges[:-1], target_edges[1:])):
        overlap = np.maximum(0.0, np.minimum(intervals[:, 1], right) - np.maximum(intervals[:, 0], left))
        if not np.any(overlap > 0):
            continue
        pose_rows = (overlap > 0) & sample.vision_observed[:, 17:20].all(axis=1)
        if pose_rows.any():
            vector = _so3_weighted_mean(raw[pose_rows, 17:20], overlap[pose_rows])
            if vector is not None:
                values[target_index, 17:20] = vector
                masks[target_index, 17:20] = True
        for start in (23, 26):
            valid_rows = (overlap > 0) & sample.vision_observed[:, start:start + 3].all(axis=1)
            if valid_rows.any():
                weights = overlap[valid_rows].astype(np.float64)
                vector = np.average(raw[valid_rows, start:start + 3], axis=0, weights=weights)
                norm = float(np.linalg.norm(vector))
                if norm > 1e-6:
                    values[target_index, start:start + 3] = vector / norm
                    masks[target_index, start:start + 3] = True
                else:
                    values[target_index, start:start + 3] = 0.0
                    masks[target_index, start:start + 3] = False
    return values, masks, coverage


def _project_posterior_text(sample: Sample, target_edges: np.ndarray) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    words, dim = sample.text_features.shape
    output = np.zeros((len(target_edges) - 1, dim), dtype=np.float32)
    mask = np.zeros((len(target_edges) - 1, dim), dtype=np.bool_)
    coverage = np.zeros((len(target_edges) - 1, dim), dtype=np.float32)
    if not words or not len(sample.ctc_times):
        return output, mask, coverage
    frame_step = CTC_FRAME_STEP_S
    first_edge = max(0.0, float(sample.ctc_times[0]) - frame_step / 2)
    frame_edges = np.clip(
        first_edge + np.arange(len(sample.ctc_times) + 1, dtype=np.float64) * frame_step,
        0.0,
        sample.duration_s,
    )
    occupancy = sample.ctc_occupancy.astype(np.float64)
    for target_index, (left, right) in enumerate(zip(target_edges[:-1], target_edges[1:])):
        overlap = np.maximum(0.0, np.minimum(frame_edges[1:], right) - np.maximum(frame_edges[:-1], left))
        if not np.any(overlap > 0):
            continue
        alpha = (occupancy * overlap[:, None]).sum(axis=0)
        alpha *= sample.text_valid
        total = float(alpha.sum())
        if total <= 0:
            continue
        output[target_index] = (alpha @ sample.text_features.astype(np.float64) / total).astype(np.float32)
        mask[target_index] = True
        coverage[target_index] = min(1.0, total / max(right - left, 1e-8))
    return output, mask, coverage


def _normalized_hard_quality(sample: Sample) -> np.ndarray:
    text_quality = np.zeros_like(sample.hard_quality, dtype=np.float32)
    valid_quality = sample.hard_valid & (sample.hard_quality > 0)
    if valid_quality.any():
        quality_reference = float(np.median(sample.hard_quality[valid_quality]))
        text_quality[valid_quality] = np.clip(
            sample.hard_quality[valid_quality].astype(np.float64) / quality_reference, 0.25, 4.0
        )
    return text_quality


def _make_views(sample: Sample) -> tuple[View, View, View]:
    edges = _grid_edges(sample.duration_s)
    hard_text_observed = np.broadcast_to(sample.text_valid[:, None], sample.text_features.shape).copy()
    text_quality = _normalized_hard_quality(sample)
    text_hard = _project_rows(
        sample.text_features.astype(np.float32), hard_text_observed,
        sample.hard_intervals, text_quality, edges,
    )
    audio = _project_rows(
        sample.audio_features.astype(np.float32), sample.audio_observed,
        _voronoi_intervals(sample.audio_times, sample.duration_s),
        np.ones(len(sample.audio_times), dtype=np.float32), edges,
    )
    vision_b0 = _project_vision_geometry(sample, edges, geometry_aware=False)
    vision_b1 = _project_vision_geometry(sample, edges, geometry_aware=True)
    speech = _project_rows(
        sample.speech_features.astype(np.float32),
        np.ones(sample.speech_features.shape, dtype=np.bool_),
        _voronoi_intervals(sample.ctc_times, sample.duration_s),
        np.ones(len(sample.ctc_times), dtype=np.float32), edges,
    )
    post_text = _project_posterior_text(sample, edges)
    view_b0 = View(
        features={"text": text_hard[0], "audio": audio[0], "vision": vision_b0[0], "speech": speech[0]},
        observed={"text": text_hard[1], "audio": audio[1], "vision": vision_b0[1], "speech": speech[1]},
        coverage={"text": text_hard[2], "audio": audio[2], "vision": vision_b0[2], "speech": speech[2]},
        edges=edges,
    )
    view_b1 = View(
        features={"text": text_hard[0], "audio": audio[0], "vision": vision_b1[0], "speech": speech[0]},
        observed={"text": text_hard[1], "audio": audio[1], "vision": vision_b1[1], "speech": speech[1]},
        coverage={"text": text_hard[2], "audio": audio[2], "vision": vision_b1[2], "speech": speech[2]},
        edges=edges,
    )
    # B2 changes only text's time projection; media timestamps and their masks stay fixed.
    view_b2 = View(
        features={**view_b1.features, "text": post_text[0]},
        observed={**view_b1.observed, "text": post_text[1]},
        coverage={**view_b1.coverage, "text": post_text[2]},
        edges=edges,
    )
    return view_b0, view_b1, view_b2


def _fit_feature_scaler(samples: list[Sample], views: list[View], indices: np.ndarray, mode: str) -> dict[str, tuple[np.ndarray, np.ndarray]]:
    modalities = ("text", "audio", "vision", "speech")
    scalers: dict[str, tuple[np.ndarray, np.ndarray]] = {}
    for modality in modalities:
        dims = views[int(indices[0])].features[modality].shape[1]
        columns: list[np.ndarray] = []
        for dim in range(dims):
            values = [
                views[int(index)].features[modality][:, dim][views[int(index)].observed[modality][:, dim]]
                for index in indices
            ]
            nonempty = [array for array in values if len(array)]
            columns.append(np.concatenate(nonempty).astype(np.float64) if nonempty else np.empty(0, dtype=np.float64))
        center = np.zeros(dims, dtype=np.float64)
        scale = np.ones(dims, dtype=np.float64)
        for dim, values in enumerate(columns):
            if not len(values):
                continue
            if mode == "robust":
                center[dim] = np.median(values)
                mad = 1.4826 * np.median(np.abs(values - center[dim]))
                std = np.std(values)
                scale[dim] = mad if mad > 1e-8 else std if std > 1e-8 else 1.0
            else:
                center[dim] = np.mean(values)
                std = np.std(values)
                scale[dim] = std if std > 1e-8 else 1.0
        scalers[modality] = (center.astype(np.float32), scale.astype(np.float32))
    return scalers


def _scaled_view(view: View, scalers: dict[str, tuple[np.ndarray, np.ndarray]], clip: bool) -> View:
    features: dict[str, np.ndarray] = {}
    for modality, value in view.features.items():
        center, scale = scalers[modality]
        transformed = (value - center[None, :]) / scale[None, :]
        if clip:
            transformed = np.clip(transformed, -8.0, 8.0)
        transformed = np.where(view.observed[modality], transformed, 0.0)
        features[modality] = transformed.astype(np.float32)
    return View(features=features, observed=view.observed, coverage=view.coverage, edges=view.edges)


def _path_signature(times: np.ndarray, values: np.ndarray) -> np.ndarray:
    # Time-augmented piecewise-linear path; return first level and antisymmetric area.
    dimension = values.shape[1] + 1
    if len(times) < 3:
        return np.zeros(dimension + dimension * (dimension - 1) // 2, dtype=np.float32)
    local_time = (times - times[0]) / max(times[-1] - times[0], 1e-8)
    path = np.column_stack((local_time, values)).astype(np.float64)
    increments = np.diff(path, axis=0)
    first = increments.sum(axis=0)
    areas: list[float] = []
    for left in range(dimension):
        for right in range(left + 1, dimension):
            prefix_left = np.cumsum(increments[:, left]) - increments[:, left]
            prefix_right = np.cumsum(increments[:, right]) - increments[:, right]
            areas.append(float(0.5 * np.sum(prefix_left * increments[:, right] - prefix_right * increments[:, left])))
    return np.concatenate((first, np.asarray(areas, dtype=np.float64))).astype(np.float32)


def _signature_for_segment(
    times: np.ndarray, values: np.ndarray, valid: np.ndarray, left: float, right: float,
) -> tuple[np.ndarray, float]:
    dim = values.shape[1]
    output_dim = (dim + 1) + (dim + 1) * dim // 2
    selected = (times >= left) & (times < right) & valid
    indices = np.flatnonzero(selected)
    if len(indices) < 3:
        return np.zeros(output_dim, dtype=np.float32), 0.0
    # Split at every missing source position; never connect across a missing interval.
    chunks = np.split(indices, np.flatnonzero(np.diff(indices) > 1) + 1)
    signatures: list[np.ndarray] = []
    weights: list[float] = []
    for chunk in chunks:
        if len(chunk) < 3:
            continue
        duration = float(max(times[chunk[-1]] - times[chunk[0]], 0.0))
        if duration <= 0:
            continue
        signatures.append(_path_signature(times[chunk], values[chunk]))
        weights.append(duration)
    if not weights:
        return np.zeros(output_dim, dtype=np.float32), 0.0
    weight_array = np.asarray(weights, dtype=np.float64)
    signature = np.average(np.stack(signatures), axis=0, weights=weight_array)
    valid_duration = min(float(np.sum(weight_array)), max(right - left, 1e-8))
    return signature.astype(np.float32), float(np.clip(valid_duration / (right - left), 0.0, 1.0))


def _dynamic_signature(sample: Sample, scaled_view: View) -> np.ndarray:
    duration = sample.duration_s
    segments = np.linspace(0.0, duration, 6)
    centers = (scaled_view.edges[:-1] + scaled_view.edges[1:]) / 2
    audio = scaled_view.features["audio"]
    audio_mask = scaled_view.observed["audio"]
    vision = scaled_view.features["vision"]
    vision_mask = scaled_view.observed["vision"]
    # log-F0, log-energy; MediaPipe jaw-open proxy and normalized mouth aperture.
    audio_values = audio[:, [67, 66]]
    audio_valid = audio_mask[:, [67, 66]].all(axis=1)
    vision_values = vision[:, [11, 31]]
    vision_valid = vision_mask[:, [11, 31]].all(axis=1)
    result: list[float] = []
    for index in range(5):
        left, right = float(segments[index]), float(segments[index + 1])
        signature, coverage = _signature_for_segment(centers, audio_values, audio_valid, left, right)
        result.extend(signature.tolist())
        result.append(coverage)
        signature, coverage = _signature_for_segment(centers, vision_values, vision_valid, left, right)
        result.extend(signature.tolist())
        result.append(coverage)
    return np.asarray(result, dtype=np.float32)


def _pool_progress(view: View, include_speech: bool, extra: np.ndarray | None = None) -> np.ndarray:
    duration = float(view.edges[-1])
    progress_edges = np.linspace(0.0, duration, 6)
    centers = (view.edges[:-1] + view.edges[1:]) / 2
    cell_lengths = np.diff(view.edges)
    names = ["text", "audio", "vision"] + (["speech"] if include_speech else [])
    pooled: list[np.ndarray] = []
    for name in names:
        values, masks, coverage = view.features[name], view.observed[name], view.coverage[name]
        rows = []
        for left, right in zip(progress_edges[:-1], progress_edges[1:]):
            select = (centers >= left) & (centers < right)
            weights = cell_lengths[select, None] * coverage[select] * masks[select]
            denom = weights.sum(axis=0)
            good = denom > 0
            out = np.zeros(values.shape[1], dtype=np.float32)
            if good.any():
                out[good] = ((values[select] * weights).sum(axis=0)[good] / denom[good]).astype(np.float32)
            rows.append(out)
        pooled.append(np.stack(rows).reshape(-1))
    if extra is not None:
        pooled.append(extra)
    return np.concatenate(pooled).astype(np.float32)


def _view_coverage(view: View, modality: str) -> float:
    mask = view.observed[modality]
    cov = view.coverage[modality]
    if not mask.size:
        return 0.0
    return float(np.mean(cov[mask])) if mask.any() else 0.0


def _metrics(y_class: np.ndarray, p_class: np.ndarray, y_value: np.ndarray, p_value: np.ndarray) -> dict[str, float]:
    pearson = float(np.corrcoef(y_value, p_value)[0, 1]) if len(y_value) > 1 and np.std(p_value) > 1e-12 else math.nan
    return {
        "accuracy": float(accuracy_score(y_class, p_class)),
        "macro_f1": float(f1_score(y_class, p_class, labels=[0, 1, 2], average="macro", zero_division=0)),
        "mae": float(mean_absolute_error(y_value, p_value)),
        "pearson": pearson,
    }


def _write_csv(path: Path, rows: list[dict[str, Any]], fields: list[str]) -> None:
    path.parent.mkdir(parents=True, exist_ok=True)
    with path.open("w", encoding="utf-8-sig", newline="") as stream:
        writer = csv.DictWriter(stream, fieldnames=fields, extrasaction="ignore")
        writer.writeheader()
        writer.writerows(rows)


def _group_bootstrap_deltas(
    samples: list[Sample], predictions: dict[str, dict[int, tuple[int, float]]], repeats: int,
) -> list[dict[str, Any]]:
    rng = np.random.default_rng(SEED)
    groups: dict[str, list[int]] = {}
    for index, sample in enumerate(samples):
        groups.setdefault(sample.video_id, []).append(index)
    group_names = sorted(groups)
    y_class = np.asarray([sample.polarity for sample in samples], dtype=np.int64)
    y_value = np.asarray([sample.sentiment for sample in samples], dtype=np.float64)
    records: list[dict[str, Any]] = []
    reference = "B1"
    for method in ("B0", "B2", "B3", "B4"):
        point_ref = np.asarray([predictions[reference][i][0] for i in range(len(samples))], dtype=np.int64)
        point_method = np.asarray([predictions[method][i][0] for i in range(len(samples))], dtype=np.int64)
        value_ref = np.asarray([predictions[reference][i][1] for i in range(len(samples))], dtype=np.float64)
        value_method = np.asarray([predictions[method][i][1] for i in range(len(samples))], dtype=np.float64)
        for metric in ("accuracy", "macro_f1", "mae", "pearson"):
            deltas: list[float] = []
            for _ in range(repeats):
                chosen = rng.choice(group_names, size=len(group_names), replace=True)
                indices = np.asarray([index for name in chosen for index in groups[str(name)]], dtype=np.int64)
                a = _metrics(y_class[indices], point_ref[indices], y_value[indices], value_ref[indices])[metric]
                b = _metrics(y_class[indices], point_method[indices], y_value[indices], value_method[indices])[metric]
                if math.isfinite(a) and math.isfinite(b):
                    deltas.append(b - a)
            ref_metric = _metrics(y_class, point_ref, y_value, value_ref)[metric]
            method_metric = _metrics(y_class, point_method, y_value, value_method)[metric]
            records.append({
                "comparison": f"{method}-B1",
                "metric": metric,
                "reference_oof": ref_metric,
                "method_oof": method_metric,
                "delta_oof": method_metric - ref_metric if math.isfinite(method_metric) and math.isfinite(ref_metric) else math.nan,
                "group_bootstrap_ci95_low": float(np.quantile(deltas, 0.025)) if deltas else math.nan,
                "group_bootstrap_ci95_high": float(np.quantile(deltas, 0.975)) if deltas else math.nan,
                "bootstrap_repeats": len(deltas),
            })
    return records


def _evaluate(samples: list[Sample], views: list[tuple[View, View, View]], bootstrap_repeats: int) -> tuple[list[dict[str, Any]], list[dict[str, Any]], list[dict[str, Any]], list[dict[str, Any]], list[dict[str, Any]]]:
    y_class = np.asarray([sample.polarity for sample in samples], dtype=np.int64)
    y_value = np.asarray([sample.sentiment for sample in samples], dtype=np.float64)
    groups = np.asarray([sample.video_id for sample in samples])
    splitter = GroupKFold(n_splits=N_FOLDS)
    folds = list(splitter.split(np.zeros(len(samples)), y_class, groups))
    methods = ("B0", "B1", "B2", "B3", "B4")
    class_oof = {method: np.full(len(samples), -1, dtype=np.int64) for method in methods}
    value_oof = {method: np.full(len(samples), np.nan, dtype=np.float64) for method in methods}
    predictions: dict[str, dict[int, tuple[int, float]]] = {method: {} for method in methods}
    fold_rows: list[dict[str, Any]] = []
    split_rows: list[dict[str, Any]] = []
    for fold_index, (train_indices, test_indices) in enumerate(folds, start=1):
        train_groups = set(groups[train_indices].tolist())
        test_groups = set(groups[test_indices].tolist())
        if train_groups & test_groups:
            raise AssertionError("video_id leakage across grouped fold")
        for index in test_indices:
            split_rows.append({
                "sample_id": samples[index].sample_id,
                "video_id": samples[index].video_id,
                "fold": fold_index,
                "split": "valid_oof",
            })

        variant_views = {
            "B0": [item[0] for item in views],
            "B1": [item[1] for item in views],
            "B2": [item[2] for item in views],
            "B3": [item[1] for item in views],
            "B4": [item[1] for item in views],
        }
        for method in methods:
            source_views = variant_views[method]
            scale_mode = "standard" if method == "B0" else "robust"
            scalers = _fit_feature_scaler(samples, source_views, train_indices, scale_mode)
            scaled = [_scaled_view(view, scalers, clip=(scale_mode == "robust")) for view in source_views]
            add_speech = method == "B4"
            x = np.stack([
                _pool_progress(
                    scaled[int(index)],
                    include_speech=add_speech,
                    extra=_dynamic_signature(samples[int(index)], scaled[int(index)]) if method == "B3" else None,
                )
                for index in train_indices
            ])
            x_valid = np.stack([
                _pool_progress(
                    scaled[int(index)],
                    include_speech=add_speech,
                    extra=_dynamic_signature(samples[int(index)], scaled[int(index)]) if method == "B3" else None,
                )
                for index in test_indices
            ])
            classifier = LogisticRegression(C=0.05, max_iter=2500, solver="lbfgs", random_state=SEED)
            classifier.fit(x, y_class[train_indices])
            predicted_class = classifier.predict(x_valid)
            regressor = Ridge(alpha=25.0, solver="lsqr")
            regressor.fit(x, y_value[train_indices])
            predicted_value = np.clip(regressor.predict(x_valid), -3.0, 3.0)
            class_oof[method][test_indices] = predicted_class
            value_oof[method][test_indices] = predicted_value
            for index, pred_class, pred_value in zip(test_indices, predicted_class, predicted_value):
                predictions[method][int(index)] = (int(pred_class), float(pred_value))
            fold_metric = _metrics(y_class[test_indices], predicted_class, y_value[test_indices], predicted_value)
            fold_rows.append({
                "method": method,
                "fold": fold_index,
                "train_samples": int(len(train_indices)),
                "valid_samples": int(len(test_indices)),
                "train_video_groups": len(train_groups),
                "valid_video_groups": len(test_groups),
                "accuracy": fold_metric["accuracy"],
                "macro_f1": fold_metric["macro_f1"],
                "mae": fold_metric["mae"],
                "pearson": fold_metric["pearson"],
                "feature_dimension": int(x.shape[1]),
            })
        print(f"[CV] fold {fold_index}/{N_FOLDS}: train={len(train_indices)} clips/{len(train_groups)} videos; valid={len(test_indices)} clips/{len(test_groups)} videos", flush=True)

    summary_rows: list[dict[str, Any]] = []
    for method in methods:
        fold_metrics = [row for row in fold_rows if row["method"] == method]
        overall = _metrics(y_class, class_oof[method], y_value, value_oof[method])
        row: dict[str, Any] = {
            "method": method,
            "sample_count": len(samples),
            "video_group_count": len(set(groups.tolist())),
            "oof_accuracy": overall["accuracy"],
            "oof_macro_f1": overall["macro_f1"],
            "oof_mae": overall["mae"],
            "oof_pearson": overall["pearson"],
        }
        for metric in ("accuracy", "macro_f1", "mae", "pearson"):
            values = np.asarray([item[metric] for item in fold_metrics if math.isfinite(item[metric])], dtype=np.float64)
            row[f"fold_{metric}_mean"] = float(values.mean()) if len(values) else math.nan
            row[f"fold_{metric}_sd"] = float(values.std(ddof=1)) if len(values) > 1 else math.nan
        summary_rows.append(row)
    prediction_rows = []
    for index, sample in enumerate(samples):
        row: dict[str, Any] = {
            "sample_id": sample.sample_id,
            "video_id": sample.video_id,
            "clip_id": sample.clip_id,
            "true_polarity": sample.polarity,
            "true_polarity_name": ("negative", "neutral", "positive")[sample.polarity],
            "true_sentiment": sample.sentiment,
        }
        for method in methods:
            row[f"{method}_predicted_polarity"] = int(class_oof[method][index])
            row[f"{method}_predicted_sentiment"] = float(value_oof[method][index])
        prediction_rows.append(row)
    delta_rows = _group_bootstrap_deltas(samples, predictions, bootstrap_repeats)
    return summary_rows, fold_rows, prediction_rows, split_rows, delta_rows


def _sample_output_rows(samples: list[Sample], views: list[tuple[View, View, View]]) -> tuple[list[dict[str, Any]], list[dict[str, Any]]]:
    sample_rows: list[dict[str, Any]] = []
    modality_rows: list[dict[str, Any]] = []
    for sample, (view0, view1, view2) in zip(samples, views):
        finite_boundaries = [item for item in sample.boundary_summary if math.isfinite(item["start_width90_s"]) and math.isfinite(item["end_width90_s"])]
        hard_text_coverage = _view_coverage(view0, "text")
        post_text_coverage = _view_coverage(view2, "text")
        sample_rows.append({
            "sample_id": sample.sample_id,
            "video_id": sample.video_id,
            "clip_id": sample.clip_id,
            "source_video_sha256": sample.video_sha256,
            "duration_s": sample.duration_s,
            "word_count": len(sample.text_words),
            "hard_aligned_word_count": int(sample.hard_valid.sum()),
            "unlocated_word_count": int((~sample.hard_valid).sum()),
            "posterior_usable_word_count": len(finite_boundaries),
            "mean_start_interval_width90_s": float(np.mean([item["start_width90_s"] for item in finite_boundaries])) if finite_boundaries else math.nan,
            "mean_end_interval_width90_s": float(np.mean([item["end_width90_s"] for item in finite_boundaries])) if finite_boundaries else math.nan,
            "text_hard_grid_coverage": hard_text_coverage,
            "text_posterior_grid_occupancy": post_text_coverage,
            "audio_grid_coverage": _view_coverage(view1, "audio"),
            "vision_grid_coverage": _view_coverage(view1, "vision"),
            "audio_native_rows": len(sample.audio_times),
            "vision_native_rows": len(sample.vision_times),
            "status": "ok",
        })
        for modality, native_length, native_dim, view, key in (
            ("text", len(sample.text_words), 768, view1, "text"),
            ("audio", len(sample.audio_times), 74, view1, "audio"),
            ("vision", len(sample.vision_times), 35, view1, "vision"),
        ):
            coverage = view.coverage[key]
            observed = view.observed[key]
            observed_duration = float((coverage * np.diff(view.edges)[:, None] * observed).sum(axis=0).mean()) if observed.size else 0.0
            modality_rows.append({
                "sample_id": sample.sample_id,
                "video_id": sample.video_id,
                "clip_id": sample.clip_id,
                "modality": modality,
                "source_duration_s": sample.duration_s,
                "observed_duration_s_mean_dimension": observed_duration,
                "native_length": native_length,
                "native_dimension": native_dim,
                "main_grid_length": len(view.edges) - 1,
                "grid_step_s": GRID_STEP_S,
                "grid_dimension": native_dim,
                "mean_grid_coverage": _view_coverage(view, key),
                "status": "ok" if native_length else "no_native_rows",
                "source_video_sha256": sample.video_sha256,
            })
    return sample_rows, modality_rows


def _create_plot(rows: list[dict[str, Any]], output_path: Path) -> None:
    labels = [row["method"] for row in rows]
    figure, axes = plt.subplots(1, 3, figsize=(14, 4.5), constrained_layout=True)
    metrics = (("oof_macro_f1", "Macro F1 (higher is better)"), ("oof_mae", "MAE (lower is better)"), ("oof_pearson", "Pearson r (higher is better)"))
    colors = ("#718096", "#2463A6", "#6A8E3A", "#C77D27", "#8357A5")
    for axis, (column, title) in zip(axes, metrics):
        values = [float(row[column]) for row in rows]
        bars = axis.bar(labels, values, color=colors)
        axis.set_title(title)
        axis.grid(axis="y", alpha=0.25)
        for bar, value in zip(bars, values):
            axis.text(bar.get_x() + bar.get_width() / 2, bar.get_height(), f"{value:.3f}", ha="center", va="bottom", fontsize=8)
    figure.suptitle("B0-B4 emotion probe on grouped out-of-fold predictions")
    output_path.parent.mkdir(parents=True, exist_ok=True)
    figure.savefig(output_path, dpi=170)
    plt.close(figure)


def _create_typical_alignment_figure(
    samples: list[Sample], views: list[tuple[View, View, View]], output_path: Path,
) -> str:
    eligible = [
        index for index, sample in enumerate(samples)
        if sample.text_words and sample.hard_valid.all() and sample.vision_observed.any()
    ]
    if not eligible:
        return ""
    median_duration = float(np.median([samples[index].duration_s for index in eligible]))
    chosen_index = min(eligible, key=lambda index: abs(samples[index].duration_s - median_duration))
    sample = samples[chosen_index]
    view_b0, _, view_b2 = views[chosen_index]
    middle = len(sample.text_words) // 2
    left_word = max(0, middle - 5)
    right_word = min(len(sample.text_words), left_word + 10)
    left_word = max(0, right_word - 10)
    word_indices = np.arange(left_word, right_word)

    figure = plt.figure(figsize=(15, 11), constrained_layout=True)
    grid = figure.add_gridspec(4, 1, height_ratios=(2.6, 1.1, 1.4, 2.1))
    alignment_axis = figure.add_subplot(grid[0, 0])
    coverage_axis = figure.add_subplot(grid[1, 0])
    feature_axis = figure.add_subplot(grid[2, 0])
    image_grid = grid[3, 0].subgridspec(1, 5, wspace=0.04)
    image_axes = [figure.add_subplot(image_grid[0, i]) for i in range(5)]

    for y_position, word_index in enumerate(word_indices):
        hard_start, hard_end = sample.hard_intervals[word_index]
        alignment_axis.plot((hard_start, hard_end), (y_position, y_position), color="#2563a6", lw=5, solid_capstyle="butt")
        boundary = sample.boundary_summary[word_index]
        if math.isfinite(boundary["start_p05_s"]) and math.isfinite(boundary["start_p95_s"]):
            center = boundary["start_mean_s"]
            alignment_axis.errorbar(
                center, y_position + 0.19,
                xerr=np.asarray([[max(0.0, center - boundary["start_p05_s"]),], [max(0.0, boundary["start_p95_s"] - center),]]),
                fmt="o", color="#d17a00", markersize=3, capsize=2, lw=1,
            )
        if math.isfinite(boundary["end_p05_s"]) and math.isfinite(boundary["end_p95_s"]):
            center = boundary["end_mean_s"]
            alignment_axis.errorbar(
                center, y_position - 0.19,
                xerr=np.asarray([[max(0.0, center - boundary["end_p05_s"]),], [max(0.0, boundary["end_p95_s"] - center),]]),
                fmt="s", color="#638b32", markersize=3, capsize=2, lw=1,
            )
    alignment_axis.set_yticks(np.arange(len(word_indices)), labels=[sample.text_words[index] for index in word_indices], fontsize=8)
    alignment_axis.set_ylim(-0.65, len(word_indices) - 0.35)
    start_values = [float(sample.hard_intervals[index, 0]) for index in word_indices]
    end_values = [float(sample.hard_intervals[index, 1]) for index in word_indices]
    alignment_axis.set_xlim(max(0.0, min(start_values) - 0.25), min(sample.duration_s, max(end_values) + 0.25))
    alignment_axis.set_xlabel("seconds from clip start")
    alignment_axis.set_title("Hard CTC word interval (blue); 90% internal start/end boundary intervals (orange/green)")
    alignment_axis.grid(axis="x", alpha=0.25)

    times = (view_b0.edges[:-1] + view_b0.edges[1:]) / 2
    coverage_axis.plot(times, view_b0.coverage["text"].mean(axis=1), label="hard word coverage", color="#2563a6")
    coverage_axis.plot(times, view_b2.coverage["text"].mean(axis=1), label="posterior word occupancy", color="#d17a00")
    coverage_axis.plot(times, view_b0.coverage["audio"].mean(axis=1), label="audio coverage", color="#638b32")
    coverage_axis.plot(times, view_b0.coverage["vision"].mean(axis=1), label="face feature coverage", color="#8357a5")
    coverage_axis.set_ylim(0, 1.05)
    coverage_axis.set_xlim(0, sample.duration_s)
    coverage_axis.set_ylabel("coverage / occupancy")
    coverage_axis.legend(ncol=4, fontsize=8, loc="upper right")
    coverage_axis.grid(alpha=0.2)

    audio = sample.audio_features.astype(np.float32)
    audio_time = sample.audio_times
    feature_axis.plot(audio_time, audio[:, 66], color="#2563a6", lw=1, label="log energy")
    pitch_valid = sample.audio_observed[:, 67]
    if pitch_valid.any():
        pitch_axis = feature_axis.twinx()
        pitch_axis.plot(audio_time[pitch_valid], np.exp(audio[pitch_valid, 67]), color="#d17a00", lw=0.8, alpha=0.8, label="F0 (Hz)")
        pitch_axis.set_ylabel("F0 (Hz)")
    vision_valid = sample.vision_observed[:, 31]
    if vision_valid.any():
        feature_axis.scatter(sample.vision_times[vision_valid], sample.vision_features[vision_valid, 31].astype(np.float32), s=9, color="#638b32", label="mouth aperture")
    feature_axis.set_xlim(0, sample.duration_s)
    feature_axis.set_xlabel("seconds from clip start")
    feature_axis.set_ylabel("log energy / mouth aperture")
    feature_axis.set_title("Example native acoustic and facial feature trajectories")
    feature_axis.grid(alpha=0.2)
    handles, labels = feature_axis.get_legend_handles_labels()
    if len(handles):
        feature_axis.legend(handles, labels, loc="upper right", fontsize=8)

    video_path = DATA_DIR / sample.video_id / f"{sample.clip_id}.mp4"
    capture = cv2.VideoCapture(str(video_path))
    frame_times = np.linspace(0.0, max(sample.duration_s - 0.05, 0.0), len(image_axes))
    for axis, timestamp in zip(image_axes, frame_times):
        capture.set(cv2.CAP_PROP_POS_MSEC, float(timestamp) * 1000.0)
        ok, bgr = capture.read()
        if ok:
            axis.imshow(cv2.cvtColor(bgr, cv2.COLOR_BGR2RGB))
            axis.set_title(f"{timestamp:.2f}s", fontsize=8)
        else:
            axis.text(0.5, 0.5, "frame unavailable", ha="center", va="center")
        axis.axis("off")
    capture.release()
    figure.suptitle(f"Typical sample alignment audit: {sample.sample_id} ({sample.duration_s:.2f}s)", fontsize=14)
    output_path.parent.mkdir(parents=True, exist_ok=True)
    figure.savefig(output_path, dpi=170)
    plt.close(figure)
    return sample.sample_id


def _write_report(
    output_dir: Path, summary: list[dict[str, Any]], deltas: list[dict[str, Any]],
    sample_rows: list[dict[str, Any]], modality_rows: list[dict[str, Any]],
    manifest: dict[str, Any],
) -> None:
    best_f1 = max(summary, key=lambda row: row["oof_macro_f1"])
    best_mae = min(summary, key=lambda row: row["oof_mae"])
    b1 = next(row for row in summary if row["method"] == "B1")
    delta = {(row["comparison"], row["metric"]): row for row in deltas}
    b4_f1 = delta[("B4-B1", "macro_f1")]
    b0_mae = delta[("B0-B1", "mae")]
    b0_pearson = delta[("B0-B1", "pearson")]
    b2_f1 = delta[("B2-B1", "macro_f1")]
    b2_mae = delta[("B2-B1", "mae")]
    b3_f1 = delta[("B3-B1", "macro_f1")]
    mean_start_width = float(np.nanmean([row["mean_start_interval_width90_s"] for row in sample_rows]))
    mean_end_width = float(np.nanmean([row["mean_end_interval_width90_s"] for row in sample_rows]))
    mean_hard_coverage = float(np.mean([row["text_hard_grid_coverage"] for row in sample_rows]))
    mean_posterior_occupancy = float(np.mean([row["text_posterior_grid_occupancy"] for row in sample_rows]))
    mean_vision_coverage = float(np.mean([row["vision_grid_coverage"] for row in sample_rows]))
    complete_hard_count = sum(row["unlocated_word_count"] == 0 for row in sample_rows)
    consistent = manifest["inputs"]["label_polarity_consistency_count"]
    lines = [
        "# 闂涓€ B0-B4 妯″瀷瀵规瘮缁撴灉",
        "",
        "鏈姤鍛婃寜 `final/Q1/闂涓€.pdf` 绗?4.8.3 鑺傚疄鏂戒簲缁勫崟鍥犵礌瀵圭収銆傚叏閮?100 鏉￠檮浠朵竴瑙嗛鍧囦娇鐢ㄥ悓涓€濂楀喕缁撶壒寰佹娊鍙栧櫒銆?.1 绉掍富鏃堕棿缃戞牸銆佽娴嬫帺鐮佷笌鎸?`video_id` 鍒嗙粍鐨勪簲鎶樺垝鍒嗐€傛儏鎰熸爣绛句粎鐢ㄤ簬鎶樺唴鐨勮交閲忛娴嬫帰閽堬紝涓嶈繘鍏ユ椂闂村榻愩€丆TC 鍚庨獙鎴栧姩鎬佺鍚嶈绠椼€?,
        "",
        "## 瀹為檯鐗瑰緛鍜屾瘮杈冪増鏈?,
        "",
        "| 鐗堟湰 | 瀹炴柦鍐呭 |",
        "| --- | --- |",
        "| B0 | CTC Viterbi 璇嶅尯闂翠笌闊抽/瑙嗚鐗╃悊鏃堕棿鍖洪棿鎶曞奖锛涜缁冩姌鍧囧€?鏍囧噯宸缉鏀撅紱瑙嗚鏃嬭浆鍚戦噺鍜岃绾垮悜閲忎綔鏅€氬垎閲忓潎鍊笺€?|",
        "| B1 | 鍦?B0 涓婃敼鐢ㄨ缁冩姌涓綅鏁颁笌 MAD 鐨勭ǔ鍋ュ昂搴︼紙MAD 涓洪浂鏃跺洖閫€鍒拌缁冩姌鏍囧噯宸級锛屽苟瀵硅瑙夋棆杞敤 SO(3) 鍧囧€笺€佽绾跨敤褰掍竴鍖栧悜閲忓潎鍊笺€?|",
        "| B2 | 鍦?B1 涓婂彧灏嗘枃鏈瘝鍚戦噺鐨勭‖杈圭晫鎶曞奖鏇挎崲涓哄浐瀹氳浆鍐?CTC 鐘舵€佸浘鐨勫墠鍚戔€撳悗鍚戝崰鎹鐜囨姇褰憋紱濯掍綋鏃堕棿鎴充粛鎸夌墿鐞嗘椂闂淬€?|",
        "| B3 | 鍦?B1 涓婁粎闄勫姞闊抽 log-F0/鑳介噺涓庤瑙?jaw-open/鍢撮儴寮€鍚堢巼鐨勬椂闂村骞夸竴銆佷簩闃惰矾寰勭鍚嶏紱缂烘祴鐐逛箣闂翠笉杩炵嚎銆?|",
        "| B4 | 鍦?B1 涓婁粎闄勫姞鍐荤粨 Wav2Vec2-base-960h 鏈€鍚庡洓灞傜殑鏃堕棿鎶曞奖琛ㄧず銆?|",
        "",
        "闄勪欢涓€鐜版湁鏍囩琛ㄥ拰鍘熷瑙嗛琚洿鎺ヤ娇鐢ㄣ€傛枃鏈负 768 缁?BERT-base-uncased 鏈洓灞傚潎鍊硷紱闊抽涓?74 缁达紙log-Mel 40銆丮FCC 13銆佄擬FCC 13銆侀煹寰?璋辩粺璁?8锛夛紱瑙嗚涓?35 缁达紙17 涓?MediaPipe blendshape 浠ｇ悊銆? 缁村ご濮裤€? 缁磋繎浼艰绾裤€? 缁撮潰閮ㄦ瘮渚嬪嚑浣曪級銆傝瑙変唬鐞嗗拰鍑犱綍绱㈠紩瀹氫箟瑙?`compare_models.py`锛屽畠浠笌 OpenFace AU 瀹氫箟骞朵笉绛夊悓锛涜闄愬埗椤诲湪璁烘枃涓槑绀恒€?,
        "",
        "## 缁撴灉",
        "",
        "鍒嗙被鎸夎繛缁爣绛剧殑涓ユ牸绗﹀彿鏋勯€?Negative/Neutral/Positive锛屽己搴﹂娴嬮檺骞呭埌 [-3, 3]銆傛瘡鎶樿缁冩姌鍗曠嫭鎷熷悎鏍囧噯鍖栧櫒銆侀€昏緫鍥炲綊锛圕=0.05锛変笌 Ridge锛坅lpha=25锛夛紱鍥哄畾浜旀姌鐢?GroupKFold 鎸夊師濮?`video_id` 鍒嗙粍銆傛€讳綋 OOF 鎸囨爣鎸?100 鏉＄暀缁勯娴嬭绠椼€?,
        "",
        "| 鏂规硶 | OOF Accuracy | OOF Macro-F1 | OOF MAE | OOF Pearson | Macro-F1 鎶樺潎鍊悸盨D |",
        "| --- | ---: | ---: | ---: | ---: | ---: |",
    ]
    for row in summary:
        lines.append(
            f"| {row['method']} | {row['oof_accuracy']:.3f} | {row['oof_macro_f1']:.3f} | {row['oof_mae']:.3f} | {row['oof_pearson']:.3f} | {row['fold_macro_f1_mean']:.3f} 卤 {row['fold_macro_f1_sd']:.3f} |"
        )
    lines += [
        "",
        f"鎸夊綋鍓?OOF 鎺㈤拡锛孧acro-F1 鏈€楂樼殑鏄?{best_f1['method']}锛坽best_f1['oof_macro_f1']:.3f}锛夛紝MAE 鏈€浣庣殑鏄?{best_mae['method']}锛坽best_mae['oof_mae']:.3f}锛夈€傜浉瀵?B1锛孊4 鐨?Macro-F1 宸负 {b4_f1['delta_oof']:+.3f}锛岃棰戠粍 Bootstrap 95% 鍖洪棿 [{b4_f1['group_bootstrap_ci95_low']:+.3f}, {b4_f1['group_bootstrap_ci95_high']:+.3f}]锛屽尯闂磋法杩囬浂銆侭0 鐨?MAE 宸负 {b0_mae['delta_oof']:+.3f}锛屽尯闂?[{b0_mae['group_bootstrap_ci95_low']:+.3f}, {b0_mae['group_bootstrap_ci95_high']:+.3f}]锛孭earson 宸负 {b0_pearson['delta_oof']:+.3f}銆侭2 鐨?Macro-F1 宸负 {b2_f1['delta_oof']:+.3f}锛屽尯闂?[{b2_f1['group_bootstrap_ci95_low']:+.3f}, {b2_f1['group_bootstrap_ci95_high']:+.3f}]锛孧AE 宸负 {b2_mae['delta_oof']:+.3f}锛汢3 鐨?Macro-F1 宸负 {b3_f1['delta_oof']:+.3f}锛屽尯闂?[{b3_f1['group_bootstrap_ci95_low']:+.3f}, {b3_f1['group_bootstrap_ci95_high']:+.3f}]銆傜粍 Bootstrap 鍙弽鏄犲綋鍓?37 涓潵婧愯棰戜笂鐨勬娊鏍蜂笉纭畾鎬э紝鏈綔澶氶噸姣旇緝鏍℃銆?,
        "",
        "## 瑕嗙洊涓庡榻愭牳楠?,
        "",
        f"鏍锋湰鏁帮細{len(sample_rows)}锛涙爣绛捐〃杩炵画鏍囩涓庢瀬鎬т竴鑷?{consistent}/100锛汣TC 纭榻愬畬鏁寸殑鏍锋湰 {complete_hard_count}/100銆侭0 纭尯闂村钩鍧囨枃鏈鐩栫巼涓?{mean_hard_coverage:.3f}锛孊2 骞冲潎璇嶅崰鎹悗楠岃川閲忎负 {mean_posterior_occupancy:.3f}锛堣繖鏄湡鏈涘崰鎹巼锛屼笉鏄鐩栫巼锛夛紱瑙嗚妫€娴嬭鐩栫巼锛堟寜 0.1 绉掔綉鏍煎钩鍧囷級涓?{mean_vision_coverage:.3f}銆傚悗楠屽唴閮?90% 杈圭晫鍖洪棿骞冲潎瀹藉害涓鸿捣鐐?{mean_start_width:.3f} 绉掋€佺粓鐐?{mean_end_width:.3f} 绉掋€傛病鏈変汉宸ヨ竟鐣屽瓙闆嗭紝鍥犳杩欎簺瀹藉害鍙槸妯″瀷鍐呴儴涓嶇‘瀹氭€ф憳瑕侊紝涓嶆槸缁忛獙鏍″噯鐜囨垨杈圭晫璇樊銆?,
        "",
        "B0鈥揃4 鐨勭壒寰佹潵婧愩€佹湁鏁堟帺鐮併€佽繛缁鐩栫巼銆佽瘝杈圭晫鍚庨獙瀹藉害銆佹姌鍐呴娴嬪強瑙嗛鍝堝笇鍧囧垎鍒繚瀛樺湪閰嶅 CSV/JSON 涓€傛病鏈変汉宸ヨ瘝杈圭晫鏃朵笉鎶ュ憡 IoU/MATE锛屼篃涓嶆妸涓嶅悓鏂规硶鐨勯娴嬫帰閽堝垎鏁拌В閲婁负瀵归綈鐪熷€笺€?,
        "",
        "## 鏂囦欢",
        "",
        "- `comparison_summary.csv`锛氫簲绉嶆柟娉曟€讳綋 OOF 涓庢姌鍧囧€兼寚鏍囥€?,
        "- `fold_metrics.csv`锛氭瘡鎶樺垎绫?鍥炲綊鎸囨爣涓庣壒寰佺淮鏁般€?,
        "- `oof_predictions.csv`锛氭瘡鏉℃牱鏈殑鐣欑粍棰勬祴銆?,
        "- `group_bootstrap_deltas.csv`锛氱浉瀵?B1 鐨勮棰戠粍閰嶅 Bootstrap 95% 鍖洪棿銆?,
        "- `sample_alignment_summary.csv`锛?00 鏉℃牱鏈€佺‖/姒傜巼鏂囨湰瑕嗙洊鍜屽悗楠岃竟鐣屽搴︺€?,
        "- `modality_summary.csv`锛?00 琛屾牱鏈?妯℃€佹槑缁嗐€?,
        "- `word_alignment_posterior.csv`锛氶€愯瘝纭竟鐣屻€佸疄闄呯浉瀵硅川閲忔潈閲嶅強鍚庨獙杈圭晫鍖洪棿銆?,
        "- `split_assignments.csv`锛氶€愭牱鏈姌鍙枫€?,
        "- `comparison.png`锛氭牳蹇?OOF 鎸囨爣鍥俱€?,
        "- `typical_alignment_example.png`锛氫腑浣嶆椂闀挎牱鏈殑璇嶈竟鐣屻€佽鐩栫巼銆佸０瀛?瑙嗚杞ㄨ抗涓庤棰戝抚鏍搁獙鍥俱€?,
        "- `run_manifest.json`锛氱幆澧冦€佹ā鍨?revision銆佸弬鏁般€佸搱甯屼笌杩愯鏃堕棿銆?,
        "",
        "## 澶嶇幇",
        "",
        "鍦ㄤ粨搴撴牴鐩綍鎵ц锛?,
        "",
        "```bash",
        "cd math",
        "uv sync",
        "uv run python compare_models.py",
        "```",
        "",
        "鐗瑰緛缂撳瓨鍐欏湪 `final/Q1/cache/native/`锛屾寮忕粨鏋滃彧鍐欏湪 `final/Q1/results/model_comparison/`銆傚垹闄ょ紦瀛樺悗浼氫粠闄勪欢涓€閲嶆柊鎻愬彇銆傞娆¤繍琛岄渶瑕佷笅杞?BERT銆乄av2Vec2 鍜?MediaPipe Face Landmarker 鏉冮噸銆?,
    ]
    (output_dir / "report.md").write_text("\n".join(lines) + "\n", encoding="utf-8")


def _write_boundary_detail(output_dir: Path, samples: list[Sample]) -> None:
    rows: list[dict[str, Any]] = []
    for sample in samples:
        quality_weights = _normalized_hard_quality(sample)
        for word_index, (word, interval, valid, quality, posterior) in enumerate(zip(
            sample.text_words,
            sample.hard_intervals,
            sample.hard_valid,
            sample.hard_quality,
            sample.boundary_summary,
        )):
            rows.append({
                "sample_id": sample.sample_id,
                "video_id": sample.video_id,
                "clip_id": sample.clip_id,
                "word_index": word_index,
                "word": word,
                "hard_start_s": float(interval[0]) if valid else math.nan,
                "hard_end_s": float(interval[1]) if valid else math.nan,
                "hard_valid": bool(valid),
                "ctc_quality_score_uncalibrated": float(quality),
                "relative_quality_weight": float(quality_weights[word_index]) if valid else math.nan,
                **posterior,
            })
    _write_csv(output_dir / "word_alignment_posterior.csv", rows, [
        "sample_id", "video_id", "clip_id", "word_index", "word", "hard_start_s", "hard_end_s",
        "hard_valid", "ctc_quality_score_uncalibrated", "relative_quality_weight", "start_mean_s", "end_mean_s",
        "start_p05_s", "start_p95_s", "end_p05_s", "end_p95_s",
        "start_width90_s", "end_width90_s",
    ])


def run(args: argparse.Namespace) -> None:
    output_dir = Path(args.output_dir).resolve()
    cache_dir = Path(args.cache_dir).resolve()
    input_dir = Path(args.data_dir).resolve()
    label_file = input_dir / "label-100.xlsx"
    records = _read_labels(label_file)
    source_hashes: dict[str, str] = {}
    for record in records:
        video_path = input_dir / record["video_id"] / f"{record['clip_id']}.mp4"
        if not video_path.is_file():
            raise FileNotFoundError(f"missing source video for {record['sample_id']}: {video_path}")
        source_hashes[record["sample_id"]] = _sha256(video_path)
    output_dir.mkdir(parents=True, exist_ok=True)
    cache_dir.mkdir(parents=True, exist_ok=True)
    cache_schema = "q1-b0b4-v3"
    face_hash = _ensure_face_model()
    samples: list[Sample] = []
    missing_records: list[dict[str, Any]] = []
    cache_paths = {
        record["sample_id"]: cache_dir / f"{_safe_name(record['video_id'])}__{_safe_name(record['clip_id'])}.npz"
        for record in records
    }
    if not args.force_extract:
        for record in records:
            path = cache_paths[record["sample_id"]]
            if path.is_file():
                try:
                    samples.append(_load_cache(path, record, source_hashes[record["sample_id"]], cache_schema))
                except (ValueError, KeyError, OSError):
                    missing_records.append(record)
            else:
                missing_records.append(record)
    else:
        missing_records = records.copy()
    # Restore source-row order even when resuming a partial cache.
    cached_by_id = {sample.sample_id: sample for sample in samples}
    if missing_records:
        device = torch.device(args.device if args.device != "auto" else ("cuda" if torch.cuda.is_available() else "cpu"))
        print(f"Extracting {len(missing_records)} uncached clips on {device}.", flush=True)
        models = _load_models(device)
        for index, record in enumerate(missing_records, start=1):
            started = time.perf_counter()
            sample = _extract_native_sample(
                record, models, device, FACE_MODEL_PATH, source_hashes[record["sample_id"]]
            )
            _save_cache(cache_paths[record["sample_id"]], sample, cache_schema)
            cached_by_id[sample.sample_id] = sample
            print(f"[features] {index}/{len(missing_records)} {sample.sample_id}: {time.perf_counter()-started:.2f}s", flush=True)
        model_info = models[-1]
    else:
        manifest_path = output_dir / "run_manifest.json"
        prior = json.loads(manifest_path.read_text(encoding="utf-8")) if manifest_path.is_file() else {}
        model_info = prior.get("models", {"text_id": TEXT_MODEL_ID, "speech_id": SPEECH_MODEL_ID, "revision_status": "features resumed from cache"})
    samples = [cached_by_id[record["sample_id"]] for record in records]
    views = [_make_views(sample) for sample in samples]
    label_consistent_count = sum(bool(record["label_consistent"]) for record in records)
    sample_rows, modality_rows = _sample_output_rows(samples, views)
    summary, folds, predictions, splits, deltas = _evaluate(samples, views, args.bootstrap_repeats)
    _write_csv(output_dir / "comparison_summary.csv", summary, [
        "method", "sample_count", "video_group_count", "oof_accuracy", "oof_macro_f1", "oof_mae", "oof_pearson",
        "fold_accuracy_mean", "fold_accuracy_sd", "fold_macro_f1_mean", "fold_macro_f1_sd",
        "fold_mae_mean", "fold_mae_sd", "fold_pearson_mean", "fold_pearson_sd",
    ])
    _write_csv(output_dir / "fold_metrics.csv", folds, [
        "method", "fold", "train_samples", "valid_samples", "train_video_groups", "valid_video_groups",
        "accuracy", "macro_f1", "mae", "pearson", "feature_dimension",
    ])
    _write_csv(output_dir / "oof_predictions.csv", predictions, [
        "sample_id", "video_id", "clip_id", "true_polarity", "true_polarity_name", "true_sentiment",
        *[f"{method}_predicted_polarity" for method in ("B0", "B1", "B2", "B3", "B4")],
        *[f"{method}_predicted_sentiment" for method in ("B0", "B1", "B2", "B3", "B4")],
    ])
    _write_csv(output_dir / "group_bootstrap_deltas.csv", deltas, [
        "comparison", "metric", "reference_oof", "method_oof", "delta_oof",
        "group_bootstrap_ci95_low", "group_bootstrap_ci95_high", "bootstrap_repeats",
    ])
    _write_csv(output_dir / "split_assignments.csv", sorted(splits, key=lambda row: row["sample_id"]), ["sample_id", "video_id", "fold", "split"])
    _write_csv(output_dir / "sample_alignment_summary.csv", sample_rows, list(sample_rows[0]))
    _write_csv(output_dir / "modality_summary.csv", modality_rows, list(modality_rows[0]))
    _write_boundary_detail(output_dir, samples)
    _create_plot(summary, output_dir / "comparison.png")
    typical_sample_id = _create_typical_alignment_figure(samples, views, output_dir / "typical_alignment_example.png")
    manifest = {
        "created_at_local": time.strftime("%Y-%m-%d %H:%M:%S %z"),
        "python": sys.version,
        "platform": platform.platform(),
        "uv_version": subprocess.run(["uv", "--version"], capture_output=True, text=True).stdout.strip(),
        "packages": {name: _version(name) for name in ("torch", "transformers", "numpy", "scipy", "scikit-learn", "mediapipe", "opencv-python-headless", "openpyxl", "matplotlib")},
        "device": args.device if args.device != "auto" else ("cuda" if torch.cuda.is_available() else "cpu"),
        "cuda_available": bool(torch.cuda.is_available()),
        "gpu_name": torch.cuda.get_device_name(0) if torch.cuda.is_available() and args.device != "cpu" else None,
        "models": model_info,
        "face_landmarker_asset_sha256": face_hash,
        "inputs": {
            "label_file": str(label_file),
            "label_file_sha256": _sha256(label_file),
            "sample_count": len(samples),
            "video_sha256_by_sample": source_hashes,
            "label_polarity_consistency_count": label_consistent_count,
        },
        "experiment": {
            "methods": {
                "B0": "quality-weighted hard CTC word projection + physical source-time pooling + train-fold standard scaling",
                "B1": "B0 + train-fold median/MAD scaling + SO(3) rotation mean + normalized gaze vector mean",
                "B2": "B1 + fixed-transcript CTC forward-backward word occupancy projection for text",
                "B3": "B1 + second-order time-augmented log signatures of selected audio/face channels",
                "B4": "B1 + frozen Wav2Vec2 last-four-layer speech representation",
            },
            "text_dimension": 768,
            "audio_dimension": 74,
            "vision_dimension": 35,
            "audio_features": list(AUDIO_NAMES),
            "vision_features": list(VISION_NAMES),
            "vision_feature_note": "17 MediaPipe blendshape proxies and approximate gaze/landmark geometry; not OpenFace AU labels",
            "grid_step_s": GRID_STEP_S,
            "audio_window_samples": FRAME_LENGTH,
            "audio_step_samples": FRAME_STEP,
            "audio_fft": N_FFT,
            "audio_mel_bands": MEL_COUNT,
            "vision_sample_rate_hz": VISION_RATE_HZ,
            "grouping": "GroupKFold by video_id",
            "fold_count": N_FOLDS,
            "probe": {"classifier": "LogisticRegression", "C": 0.05, "regressor": "Ridge", "alpha": 25.0, "intensity_clip": [-3.0, 3.0]},
            "bootstrap": {"unit": "video_id", "repeats": args.bootstrap_repeats, "seed": SEED},
            "human_boundary_reference_count": 0,
            "boundary_interval_level": 0.90,
            "hard_text_quality_weight": "per-sample median-normalized CTC path score clipped to [0.25, 4.0]; uncalibrated",
            "cache_schema": cache_schema,
        },
        "output_dir": str(output_dir),
        "typical_sample_id": typical_sample_id,
        "elapsed_seconds": time.time() - args.started_at,
    }
    (output_dir / "run_manifest.json").write_text(json.dumps(manifest, ensure_ascii=False, indent=2) + "\n", encoding="utf-8")
    _write_report(output_dir, summary, deltas, sample_rows, modality_rows, manifest)
    print(f"Saved B0-B4 comparison to {output_dir}", flush=True)


def _parse_args() -> argparse.Namespace:
    parser = argparse.ArgumentParser(description="Run Q1 B0-B4 model comparison on Attachment 1.")
    parser.add_argument("--data-dir", type=Path, default=DATA_DIR)
    parser.add_argument("--output-dir", type=Path, default=OUTPUT_DIR)
    parser.add_argument("--cache-dir", type=Path, default=CACHE_DIR)
    parser.add_argument("--device", choices=("auto", "cuda", "cpu"), default="auto")
    parser.add_argument("--bootstrap-repeats", type=int, default=2000)
    parser.add_argument("--force-extract", action="store_true", help="recompute native features from the 100 raw videos")
    args = parser.parse_args()
    if args.bootstrap_repeats < 100:
        parser.error("--bootstrap-repeats must be at least 100")
    args.started_at = time.time()
    return args


if __name__ == "__main__":
    run(_parse_args())