// Block bootstrap for dependent data (registry keys bootstrap.*):
//   stationary (Politis & Romano 1994): blocks start at uniform random indices, lengths ~ Geometric(p = 1/mean_block),
//                                       wrapping circularly; expected block length = mean_block.
//   circular block (Politis & Romano 1992): fixed block length b, start indices uniform on [0, n), wrap around.
//   iid (for comparison only; not valid for dependent data — flagged in the registry).
// Index generation is deterministic given (seed, replica id): replica r uses Rng(seed, r), so replicas can be generated
// in parallel (CPU threads) or on the GPU with the same streams.
#pragma once
#include <vector>
#include <span>
#include <functional>
#include "raid/core/types.hpp"
#include "raid/core/rng.hpp"
#include "raid/core/plan.hpp"
#include "raid/core/estimate.hpp"

namespace raid::stats {

enum class BootstrapScheme { stationary, circular_block, iid };

struct BootstrapOptions {
  BootstrapScheme scheme{BootstrapScheme::stationary};
  f64 mean_block{0};          // 0 → n^{1/3} (rule of thumb; the paper reports the value used)
  usize replicas{999};
  u64 seed{20261102};
  unsigned threads{0};
  f64 ci_level{0.95};
};

/// Index vector of one replica (length n).
std::vector<u32> bootstrap_indices(usize n, const BootstrapOptions& o, usize replica);

struct BootstrapResult {
  std::vector<f64> statistics;   // one per replica (NaN if the statistic failed)
  f64 point{0};                  // statistic on the original sample
  f64 se{0};                     // bootstrap standard error
  f64 ci_lo{0}, ci_hi{0};        // percentile CI
  f64 bias{0};                   // mean(stat*) − point
  usize failures{0};
  f64 block_length{0};
  f64 seconds{0};
};

/// Generic scalar statistic bootstrap: stat(x_resampled) evaluated for each replica.
BootstrapResult bootstrap(std::span<const f64> x, const std::function<f64(std::span<const f64>)>& stat, const BootstrapOptions& o);

/// Multivariate: statistic returns a vector (e.g. regression coefficients); rows of (y, X) are resampled jointly.
struct BootstrapVecResult { std::vector<std::vector<f64>> statistics; std::vector<f64> point, se, ci_lo, ci_hi; usize failures{0}; f64 block_length{0}; f64 seconds{0}; };
BootstrapVecResult bootstrap_rows(usize n, const std::function<std::vector<f64>(std::span<const u32>)>& stat_of_indices, const BootstrapOptions& o);

/// Realized-variance bootstrap: for each replica, resample intraday returns within a session by blocks and recompute
/// RV (Gonçalves & Meddahi 2009 use the iid wild bootstrap for RV; we use block schemes to respect dependence).
BootstrapResult bootstrap_rv(std::span<const f64> intraday_returns, const BootstrapOptions& o);

}  // namespace raid::stats