Небинированный событийный анализ (unbinned):
зачем работать с каждым событием

1.Зачем нужен небинированный анализ

Стандартный подход в физике высоких энергий — биновый анализ: события группируются в бины гистограммы, и для каждого бина строится пуассоновское правдоподобие. Модификаторы (NormSys, HistoSys, ShapeSys и др.) описывают систематические неопределенности. Этот подход доминирует десятилетия по веским причинам: он хорошо масштабируется, хорошо понятен и имеет зрелые реализации (ROOT/RooFit, pyhf, NextStat).

Но биннинг — это выбор, и он имеет цену. Когда вы сбрасываете события в бины, вы теряете информацию о точном значении наблюдаемой внутри бина. В пределе бесконечной статистики и бесконечно мелких бинов потеря стремится к нулю. На практике — нет.

Когда биннинг вредит

• ›Узкие резонансы. Сигнал с шириной σ ≈ 2 ГэВ при размере бина 5 ГэВ размывается. Чувствительность деградирует, потому что события из ядра пика и из хвоста оказываются в одном бине с одинаковым весом.
• ›Малая статистика. При 50–200 событиях в сигнальной области любой выбор биннинга — артефакт. Бины с 0–2 событиями плохо описываются пуассоновским приближением. Небинированное правдоподобие не вводит дополнительной дискретизации из-за выбора бинов.
• ›Многомерные наблюдаемые. При биннинге в N измерениях число бинов растёт как O(k^N). Для 3D-наблюдаемой с 50 бинами на ось потребуется 125 000 бинов — большинство пустых. Небинированный подход не вынуждает платить цену «комбинаторной пустоты» ради строгости.
• ›Shape-систематики на уровне событий. Когда вес каждого события зависит от мешающего параметра (например, PDF-неопределенности в p-p столкновениях), побинная аппроксимация огрубляет зависимость. Событийная модель учитывает каждый вес (или сдвиг наблюдаемой) явно.

Почему небинированное правдоподобие тяжёлое вычислительно

Информационная строгость unbinned-подхода оплачивается вычислениями. В биновом анализе одна оценка NLL — это цикл по бинам. В небинированном — цикл по событиям: для каждого события нужно вычислить смесь процессов Σ_p ν_p · f_p(x_i), а затем суммировать log.

В первом приближении стоимость одной оценки NLL масштабируется как O(N_events × N_processes). Для градиента добавляется стоимость производных PDF по параметрам формы. Дальше это умножается на число итераций оптимизатора и на число оценок функции внутри line search.

Практический вывод: при миллионах событий наивная реализация быстро становится непригодной. Поэтому NextStat уделяет внимание раскладке памяти (SoA), batch-оценке PDF и, для части моделей, GPU-ускорению.

2.Математические основы

2.1 Биновое правдоподобие (HistFactory)

В стандартном биновом HistFactory-анализе правдоподобие строится по бинам гистограммы:

L_binned = ∏_bins Poisson(n_b | ν_b(θ)) × ∏_NP Constraint(θ_j)

где ν_b(θ) = Σ_samples  μ_s · f_s(θ) · h_s,b
    n_b    — наблюдённое число событий в бине b
    θ      — мешающие параметры (систематики)

2.2 Небинированное расширенное правдоподобие

В небинированном подходе правдоподобие строится по каждому событию индивидуально. Для канала с несколькими процессами (сигнал + фоны) используется расширенное правдоподобие (extended likelihood):

L_unbinned = e^(-ν_total) × ∏_events [ Σ_processes  ν_p · f_p(x_i | θ) ]
           × ∏_NP Constraint(θ_j)

где ν_total = Σ_processes ν_p        — суммарное ожидаемое число событий
    ν_p                               — ожидаемое число событий процесса p
    f_p(x | θ)                        — нормализованная PDF процесса p
    x_i                               — наблюдаемые i-го события

Множитель e^(-ν_total) — расширенный (extended) терм. Он делает правдоподобие чувствительным не только к форме распределения, но и к числу событий, что необходимо для оценки параметра силы сигнала μ.

2.3 Отрицательное лог-правдоподобие (NLL)

На практике минимизируется NLL:

−log L = ν_total − Σ_events log( Σ_processes ν_p · f_p(x_i | θ) )
       + Σ_NP penalty(θ_j)

Это в точности та функция, которую NextStat передаёт в L-BFGS-B оптимизатор через trait LogDensityModel. Тот же оптимизатор, тот же профильный скан, та же проверка гипотез — как для бинов, так и для небинированных моделей.

3.Каталог функций плотности вероятности

NextStat предоставляет набор параметрических и непараметрических PDF, каждая из которых реализует trait UnbinnedPdf с двумя ключевыми методами:

• ›log_prob_batch(events, params, out) — вычисление log f(x|θ) для всех событий за один вызов
• ›log_prob_grad_batch(events, params, out_logp, out_grad) — то же + аналитический градиент ∂log f / ∂θ

Batch-интерфейс критичен для производительности: events хранятся в колоночном формате (SoA), и PDF может использовать SIMD-инструкции для параллельного вычисления по всем событиям.

3.1 Параметрические PDF для сигнала

• ›Gaussian [μ, σ] — стандартное гауссово распределение. Годится для чистых резонансов с симметричным детекторным разрешением.
• ›Crystal Ball [μ, σ, α, n] — гауссово ядро со степенным хвостом с одной стороны. Описывает детекторные эффекты (утечку энергии, мёртвый материал) в калориметрических измерениях. Стандартная PDF для H→γγ, H→ZZ*→4ℓ, J/ψ→μμ.
• ›Double Crystal Ball [μ, σ, α_l, n_l, α_r, n_r] — степенные хвосты с обеих сторон. Используется когда и left- и right-tail имеют непуассоновскую природу.

3.2 PDF для фона

• ›Exponential [λ] — e^(λx) на конечном интервале. Классическое описание комбинаторного фона.
• ›Chebyshev [c₁, …, c_k] — полином Чебышёва порядка k на [-1, 1], отображённый на поддержку наблюдаемой. Гладкая и стабильная аппроксимация фона произвольной формы. Порядок определяется числом параметров.

3.3 Шаблонные и непараметрические PDF

• ›Histogram — кусочно-постоянная PDF по бинам гистограммы. Позволяет использовать MC-шаблон как PDF в небинированном анализе: событиям присваивается плотность бина, в который они попадают.
• ›Histogram из TTree — то же, но гистограмма строится на лету из ROOT TTree.
• ›KDE (Kernel Density Estimation) — ядерная оценка плотности по набору тренировочных событий. Непараметрическая — не требует выбора аналитической формы. Адаптивная bandwidth.
• ›Morphing Histogram — гистограмма с HistoSys-систематикой формы (up/down шаблоны), интерполируемая по мешающим параметрам.
• ›Morphing KDE — KDE со событийными весовыми систематиками (weight_systematics). Позволяет описывать зависимость формы распределения от мешающих параметров через up/down веса для каждого события.
• ›KDE с горизонтальным морфингом (HorizontalMorphingKdePdf) — KDE, в которой мешающий параметр может сдвигать центры ядер (горизонтальные систематики, horizontal_systematics). Реализует нормировку на ограниченной поддержке наблюдаемой, аналитические градиенты и sample()для генерации тоев.

4.EventStore: колоночное хранилище

Небинированное правдоподобие требует вычисления f(x_i) для каждого события. При миллионах событий и десятках PDF-вызовов на каждой итерации оптимизатора, раскладка памяти критична.

EventStore использует Structure-of-Arrays (SoA): каждая наблюдаемая хранится как непрерывный вектор f64. Это обеспечивает:

• ›локальность кеша — при вычислении PDF по одной наблюдаемой данные лежат непрерывно в памяти
• ›SIMD-дружелюбность — 4 события за раз (f64x4) без перепаковки
• ›Раздельный доступ — PDF «mass» не трогает столбец «pt»

Источники данных:

{
  "data": {
    "file": "data.root",
    "tree": "CollectionTree",
    "selection": "mass > 100 && mass < 150 && nJets >= 2",
    "weight": "mcWeight * puWeight"
  },
  "observables": [
    { "name": "mass", "bounds": [100, 150] },
    { "name": "BDT_score", "expr": "bdt_out", "bounds": [0, 1] }
  ]
}

NextStat нативно читает ROOT TTree (через ns-root), поддерживает expression-based observables (произвольные выражения над ветками) и selection cuts. Никакой зависимости от ROOT C++ — только чистый Rust.

Второй поддерживаемый формат — Parquet. Это важно для репликабельности: файл можно читать и валидировать вне ROOT-экосистемы, а контракт «события → столбцы» становится прозрачным.

Минимальный контракт Parquet для событийного ввода:

- Один Parquet-файл = таблица событий (по строке на событие)
- Каждая наблюдаемая = отдельная колонка Float64 с именем, совпадающим с observable.name
- Опционально:
  - _weight (Float64): частотный вес события (по умолчанию 1.0)
  - _channel (Utf8): метка канала для multi-channel файлов
- В footer metadata (key-value) NextStat пишет:
  - nextstat.schema_version = nextstat_unbinned_events_v1
  - nextstat.observables = JSON списка наблюдаемых и их bounds

Критический инвариант: веса наблюдённых данных трактуются как неотрицательные частотные веса. На входе проверяется, что каждый w_i конечен, что w_i ≥ 0, и что Σ w_i > 0. Отрицательные веса не допускаются как контракт модели угроз: они ломают вероятностную интерпретацию расширенного правдоподобия и часто приводят к дивергенции оптимизатора.

5.UnbinnedModel: сборка правдоподобия

UnbinnedModel собирает правдоподобие из каналов, процессов и ограничений. Структура намеренно повторяет архитектуру HistFactory:

Model
├── Parameters (name, init, bounds, constraint?)
├── Channel "SR"
│   ├── data: ROOT TTree → EventStore
│   ├── observable: mass [100, 150]
│   ├── Process "signal"
│   │   ├── pdf: crystal_ball(mass_mu, mass_sig, alpha, n)
│   │   └── yield: 50 × μ, modified by NormSys(alpha_lumi)
│   └── Process "background"
│       ├── pdf: exponential(lambda)
│       └── yield: fixed(1000)
└── Channel "CR" (optional, control region)

Выражения для выхода (yield)

Каждый процесс определяет, сколько событий ожидается:

• ›Fixed(ν) — фиксированное число. Для хорошо известного фона.
• ›Parameter(index) — свободный параметр, оптимизируемый напрямую. Для data-driven оценки фона.
• ›Scaled(base, μ) — base_yield × μ. Классический signal strength.
• ›Modified(base, modifiers) — yield с мультипликативными модификаторами: NormSys и WeightSys.

Модификаторы выхода (rate modifiers)

Модификаторы yield повторяют семантику HistFactory:

• ›NormSys(α, lo, hi) — кусочно-экспоненциальная интерполяция. При α=0 фактор равен 1; при α=+1 — hi; при α=-1 — lo. Описывает нормировочные систематики (luminosity, cross-section).
• ›WeightSys(α, lo, hi, interp_code) — HistoSys-стиль интерполяция (code0 — кусочно-линейная, code4p — кусочно-экспоненциальная). Для rate-only систематик с нелинейным поведением.

6.Типичный workflow: поиск резонанса

Рассмотрим поиск узкого резонанса в спектре инвариантной массы. Сигнал описывается Crystal Ball, фон — экспонентой. 50 ожидаемых сигнальных событий, 1000 фоновых, нормировочная систематика 5%.

Шаг 1: JSON-спецификация

{
  "schema_version": "nextstat_unbinned_spec_v0",
  "model": {
    "poi": "mu",
    "parameters": [
      { "name": "mu",         "init": 1.0,   "bounds": [0, 10] },
      { "name": "mass_mu",    "init": 125.0,  "bounds": [120, 130] },
      { "name": "mass_sigma", "init": 2.0,    "bounds": [0.5, 5] },
      { "name": "alpha_cb",   "init": 1.5,    "bounds": [0.5, 5] },
      { "name": "n_cb",       "init": 3.0,    "bounds": [1, 20] },
      { "name": "lambda",     "init": -0.02,  "bounds": [-1, 0] },
      { "name": "alpha_lumi", "init": 0,      "bounds": [-5, 5],
        "constraint": { "type": "gaussian", "mean": 0, "sigma": 1 } }
    ]
  },
  "channels": [{
    "name": "SR",
    "data": { "file": "data.root", "tree": "events",
              "selection": "mass > 100 && mass < 150" },
    "observables": [{ "name": "mass", "bounds": [100, 150] }],
    "processes": [
      {
        "name": "signal",
        "pdf": { "type": "crystal_ball", "observable": "mass",
                 "params": ["mass_mu", "mass_sigma", "alpha_cb", "n_cb"] },
        "yield": { "type": "scaled", "base_yield": 50, "scale": "mu",
                   "modifiers": [{ "type": "normsys", "param": "alpha_lumi",
                                   "lo": 0.95, "hi": 1.05 }] }
      },
      {
        "name": "background",
        "pdf": { "type": "exponential", "observable": "mass",
                 "params": ["lambda"] },
        "yield": { "type": "fixed", "value": 1000 }
      }
    ]
  }]
}

Спецификация — это контракт: однозначное описание модели, данных и ограничений, которое одинаково интерпретируется локально, в CI и при репликации третьими сторонами. Файл читается как YAML по умолчанию; если расширение .json, ожидается JSON.

Чтобы избежать двусмысленностей, важно понимать, как устроены ключевые секции.

• ›schema_version. Должен быть ровно nextstat_unbinned_spec_v0. Это якорь для валидации и обратной совместимости.
• ›model.poi. Имя параметра интереса (POI). Оно необходимо для команд профилирования и тоев (например,unbinned-scan и unbinned-fit-toys): именно этот параметр сканируется и суммаризируется.
• ›model.parameters. Единый список параметров. Имена используются далее как ссылки в pdf.params, в модификаторах и в разделе yield. bounds задают область допустимых значений; init — точка старта оптимизатора.constraint добавляет штраф в NLL (например, гауссово ограничение для систематики светимости).
• ›channels. Канал — это набор данных + набор процессов, которые объясняют эти данные. У каждого канала:
  • ›data.file и data.tree указывают ROOT-файл и TTree. Относительные пути разрешаются относительно каталога конфигурационного файла.
  • ›data.selection — выражение отбора событий (cuts). Оно вычисляется на номинальных ветках; для «горизонтальных» систематик (см. ниже) отбор не пересчитывается.
  • ›observables задают имена наблюдаемых и их поддержкуbounds. Эти границы — часть контракта нормировки.
• ›processes. Процесс состоит из PDF и ожидаемого выхода (yield).
  • ›pdf задаёт тип (например,crystal_ball) и список параметров формы в params.
  • ›yield задаёт ожидаемое число событий:fixed (фиксированное),parameter (свободный параметр) илиscaled (база × POI). Модификатор normsys реализует HistFactory-подобную систематику нормы (фактор lo/hi управляется nuisance-параметром).

Для непараметрических форм, построенных прямо из дерева (histogram_from_tree и kde_from_tree), поддерживаются два типа «событийных» систематик:

• ›weight_systematics. up/down выражения задают отношения весов w_up / w_nom и w_down / w_nom на уровне событий. Они могут влиять на форму (пересборка шаблона/весов KDE) и/или на выход (yield) — в зависимости от флагов применения.
• ›horizontal_systematics. up/down выражения задают значение наблюдаемой при α=±1. Это «горизонтальный» морфинг: меняется положение события по оси наблюдаемой, но номинальные веса остаются теми же. Один nuisance-параметр не может одновременно быть источником shape через веса и через горизонтальный сдвиг.

Событийные систематики: веса и горизонтальные сдвиги

Практически полезно мыслить о двух «ручках» как о принципиально разных классах эффектов.

• ›Весовые систематики (weight_systematics) меняют вклад события в оценку плотности/шаблона.
• ›Горизонтальные систематики (horizontal_systematics) меняют «координату события» по оси наблюдаемой. Для kde_from_tree это означает морфинг центров ядер, а для histogram_from_tree — пересчёт того, как события заполняют бины.

Минимальный пример горизонтальной систематики в kde_from_tree (выражения up/down задают значение наблюдаемой при α=±1):

{
  "pdf": {
    "type": "kde_from_tree",
    "observable": "mass",
    "source": {
      "file": "mc.root",
      "tree": "events",
      "selection": "mass > 100 && mass < 150",
      "weight": "w"
    },
    "horizontal_systematics": [{
      "param": "alpha_mass_scale",
      "up": "mass * 1.001",
      "down": "mass * 0.999",
      "interp": "code4p"
    }]
  }
}

Ключевой инвариант для предсказуемости: один и тот же nuisance-параметр не должен одновременно задавать shape через weight_systematics и через horizontal_systematics. Это защищает от двойного учёта одного и того же эффекта и от неоднозначности интерпретации.

Если вы хотите проверить файл до запуска фита, распечатайте JSON-схему и провалидируйте спецификацию:

nextstat config schema --name unbinned_spec_v0

Шаг 2: MLE-фит

nextstat unbinned-fit --config search.json --output result.json

Результат содержит best-fit значения всех параметров (включая μ̂), неопределенности по гессиану, NLL в минимуме и статус сходимости.

JSON-вывод MLE-фита — это тоже контракт. Минимально он включает:input_schema_version=nextstat_unbinned_spec_v0, массивы parameter_names / bestfit /uncertainties, скаляр nll, флагconverged и счётчик итераций n_iter. Это позволяет строить автоматические отчёты и валидаторы без парсинга stdout-человеческого формата.

Корректность как гейт: closure/coverage и RooFit reference

Для unbinned-пайплайна скорость не имеет смысла без двух уровней доказательств: (1) закрытие (closure) на синтетике и (2) sanity-check против независимой реализации (RooFit) на канонических 1D моделях.

• ›Closure (быстро, по умолчанию в тестах): генерируем большой датасет из truth и проверяем, что фит восстанавливает параметры в пределахUNBINNED_CLOSURE_PARAM_ATOL=0.15 илиUNBINNED_CLOSURE_PARAM_RTOL=0.10.
• ›Coverage (медленно, opt-in): прогоняем много тоев и проверяем долю попаданий truth POI в 1σ интервал. Контракт допускает пере-покрытие (консервативный гессиан), поэтому диапазон широк: [0.55, 0.98].
• ›RooFit reference: отдельный скрипт строит тот же датасет и сравнивает NextStat с RooFit. Порог по NLL намеренно мягкий (NLL_ATOL=0.5), потому что у тулов могут отличаться константы нормировки; по параметрам используется PARAM_ATOL=0.3 / PARAM_RTOL=0.10.

# Closure (параметры должны восстановиться в пределах допусков)
pytest -v tests/python/test_unbinned_closure_coverage.py -k closure

# Coverage (slow, opt-in)
NS_RUN_SLOW=1 NS_TOYS=200 pytest -v -m slow tests/python/test_unbinned_closure_coverage.py -k coverage

# RooFit reference (требует ROOT; иначе сравнение пропускается)
python tests/validate_roofit_unbinned.py --cases gauss_exp,cb_exp --seed 42 --output tmp/validate_roofit_unbinned.json

GPU-режим (консервативный subset)

Для ускорения вычислений NextStat поддерживает GPU-оценку NLL+градиента для ограниченного подмножества unbinned-моделей. GPU включается явно флагом --gpu:

nextstat unbinned-fit --config search.json --output result.json --gpu metal
nextstat unbinned-fit --config search.json --output result.json --gpu cuda

• ›Что ускоряется. GPU берёт на себя NLL+градиент (редукции по событиям), а оптимизация остаётся тем же L-BFGS-B с теми же границами.
• ›PDF на GPU. Поддерживаются gaussian, exponential,crystal_ball, double_crystal_ball,chebyshev (order ≤ 16) и histogram. CPU-only типы (например, argus, voigtian,kde, spline) автоматически требуют CPU-путь.
• ›Yield и модификаторы. Поддерживаются fixed, parameter,scaled, а также модификаторы нормыnormsys и weightsys.
• ›Данные. GPU-путь рассчитан на 1D наблюдаемую на канал; multi-channel поддерживается суммированием вкладов каналов. Частотные веса событий поддерживаются, и проходят ту же валидацию неотрицательности, что и CPU.
• ›Точность. cuda — f64 (паритет с CPU до редукционных эффектов);metal — f32 (на Apple Silicon), поэтому допуски шире.

Если спецификация выходит за эти рамки (например, KDE/spline/voigtian, нейронные PDF без GPU-пути, или нестандартные расширения), используйте CPU-путь без --gpu. Это осознанная стратегия: сначала — узкий, проверяемый GPU subset, затем расширение покрытия.

Для toy-прогонов существует отдельный контракт: unbinned-fit-toys может использовать --gpu cuda|metal для batch-фитов, а опция--gpu-sample-toys (при наличии поддерживаемого подмножества PDF) держит toy-данные device-resident, уменьшая лишние копирования.

Шаг 3: Профильный скан μ

nextstat unbinned-scan --config search.json \
  --start 0 --stop 5 --points 51

Для каждого значения μ фиксируется POI, а все остальные параметры профилируются. Результат — кривая Δ(−2 log L) vs μ, из которой извлекается доверительный интервал.

Шаг 4: тои (coverage, bias, pulls)

nextstat unbinned-fit-toys --config search.json \
  --n-toys 1000 --seed 42 --gen init

Для каждого тоя: Poisson(ν_total) → сэмплирование из PDF → полный MLE-фит. Это позволяет проверить coverage, bias и pull distribution — критичные метрики для валидации статистической процедуры.

Для автоматических прогонов (CI) доступны гейты качества: например, требование сходимости всех тоев и ограничения на среднее/σ pull для POI (см. флаги --require-all-converged,--max-abs-poi-pull-mean, --poi-pull-std-range). При провале гейтов артефакты всё равно записываются, а команда завершается с ненулевым кодом.

Шаг 5: проверка гипотез (toy-based CLs, qtilde)

Для верхних пределов и CLs-логики в HEP обычно используют test statistic qtilde (q̃μ): это one-sided версия профильной статистики, которая обнуляется, если наблюдённая оценка μ̂ оказывается выше тестируемого μ.

В NextStat toy-based CLs для небинированных моделей реализован через sampler hookв inference-слое: алгоритм (qtilde + подсчёт хвостовых вероятностей) общий, а генерация тоев передаётся как функция sample_toy(model, params, seed). Для событийных моделей эта функция обычно сводится к model.sample_poisson_toy(params, seed).

CLI-обёртка для production-проверок:

nextstat unbinned-hypotest-toys --config search.json --mu 1.0 \
  --n-toys 1000 --seed 42 --expected-set --threads 1

Команда генерирует две toy-ансамбли (b-only при μ=0 и s+b при μ=μ_test), строит распределения qtilde и возвращает CLs, CLsb и CLb вместе с наблюдённым q_obs и μ̂. Опция --expected-set дополнительно отдаёт «бразильский набор» (5 точек expected CLs в стандартном порядке [2, 1, 0, -1, -2]).

Для корректности toy-based CLs требуется, чтобы μ=0 лежал внутри границ POI (иначе b-only гипотеза не определена). Для репликабельности генерация использует раздельные детерминированные seed’ы для b-only и s+b ансамблей.

7.Биновый и небинированный подходы: когда что использовать

В NextStat оба подхода сосуществуют и используют один оптимизатор. Переключение между ними — вопрос формата входных данных, а не смены инструмента.

8.Текущий статус и что впереди

Небинированный анализ в NextStat находится на этапе 1: полностью рабочая реализация с 10 PDF, расширенным правдоподобием (extended likelihood), CLI-командами и JSON-спецификацией.

• ›Этап 1 (сейчас) — 1D PDF, MLE-фит, профильный скан, тои, спецификация v0, консервативный GPU subset для NLL+градиента
• ›Этап 2 (следующий) — многомерные PDF (2D Gaussian, conditional PDF), расширение GPU-покрытия, Python API
• ›Этап 3 — полная совместимость с HistFactory (небинированные каналы в том же workspace, что и биновая часть)

Материал этой статьи соответствует NextStat 0.9.0.

Связанное чтение

• ›Документация: Небинированный анализ — API, каталог PDF, JSON-спецификация
• ›Документация: Модели HistFactory — биновый анализ для сравнения
• ›Справка CLI — все команды unbinned-fit / unbinned-scan / unbinned-fit-toys
• ›Как NextStat делает HistFactory дифференцируемым