commit b0a7a7cf7abdbae101d2efb87bcaf13a132b7501
parent 8c7d3c072dba28c933e84c15ddff14a8ded14509
Author: Jared Tobin <jared@jtobin.io>
Date: Sun, 28 Jun 2026 20:57:37 -0230
tighten correctness and validation across the API
Acts on a code review of v0.1.0. The headline change is that
'decide' now reflects "ever crossed the threshold," not "currently
above the threshold" -- this is the supremum-style event Ville's
inequality actually bounds. State carries a running max log-wealth
per direction (Bounded) or globally (Bernoulli), and 'log_wealth'
returns the supremum so far.
Other correctness fixes:
* 'Fixed' bettor is now clipped to '[0, lambda_max]' like the
other bettors. The Bounded module no longer needs the 'tiny'
wealth-factor floor; that was hiding a violation rather than
fixing one.
* The Newton bettor uses the WSR (2024) Algorithm 2 learning
rate '2 / (2 - log 3)', matching the paper.
API hardening:
* Config constructors return 'Either ConfigError Config' and
reject invalid statistical parameters (alpha out of (0,1),
lo >= hi, m outside (lo,hi), p0 outside (0,1)).
* 'Numeric.Eproc.Bernoulli.config' argument order changed from
(alpha, p0, bettor) to (p0, alpha, bettor) so alpha is
second-to-last across all three test modules.
Internal cleanup:
* Bettor logic (BetState, init_bet, bet_lambda, step_bet) lives
in 'Numeric.Eproc.Common'; Bounded and Bernoulli are thin
consumers. Bench numbers unchanged.
Docs:
* Module headers reframe H_0 in conditional-expectation form,
which is what anytime validity actually requires.
* 'update' docstrings state the support precondition explicitly.
Tests: 16 -> 41. Added latched-rejection cases, config-validation
group, and QuickCheck properties for log-wealth finiteness,
monotonicity, Fixed-bettor safety under arbitrary lambdas, and
the latched-rejection invariant. Tightened null-FPR slack
0.10 -> 0.08. Added a Bernoulli Fixed-bettor smoke test.
Benches: added Bernoulli update/fold groups.
Diffstat:
8 files changed, 656 insertions(+), 360 deletions(-)
diff --git a/README.md b/README.md
@@ -16,24 +16,26 @@ A sample GHCi session:
> import qualified Numeric.Eproc.Bounded as Bounded
>
> -- hypothesis: E[X] = 0.5 for samples in [0, 1] at alpha = 1e-3, tested
- > -- with the Newton bettor
- > let cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
- > let s0 = Bounded.initial cfg
+ > -- with the Newton bettor. 'config' returns 'Either ConfigError Config'
+ > -- and refuses inputs outside the mathematical regime.
+ > let Right cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ > let s0 = Bounded.initial cfg
>
> -- ten observations (drifting from hypothesis), and state afterwards
> let xs = [1, 1, 0, 1, 1, 0, 1, 1, 1, 1]
> let s10 = foldl' (Bounded.update cfg) s0 xs
>
- > -- inspect wealth and stopping decision at any point
+ > -- inspect (supremum-so-far) log-wealth and stopping decision at any
+ > -- point
> Bounded.log_wealth s10
- 0.7182493502552663
+ 0.4054651081081644
> Bounded.decide cfg s10
Continue
>
> -- with enough evidence, the hypothesis is rejected
> let s300 = foldl' (Bounded.update cfg) s0 (concat (replicate 30 xs))
> Bounded.log_wealth s300
- 53.092214534054165
+ 51.142711428622924
> Bounded.decide cfg s300
Reject
```
diff --git a/bench/Main.hs b/bench/Main.hs
@@ -4,6 +4,7 @@
module Main where
import Control.DeepSeq
+import qualified Numeric.Eproc.Bernoulli as Bern
import qualified Numeric.Eproc.Bounded as Bounded
import qualified Numeric.Eproc.Paired as P
import Criterion.Main
@@ -12,35 +13,44 @@ import Criterion.Main
-- WHNF == NF for these types. orphan instances keep the library API
-- untouched.
instance NFData Bounded.State where rnf !_ = ()
-instance NFData P.State where rnf !_ = ()
+instance NFData P.State where rnf !_ = ()
+instance NFData Bern.State where rnf !_ = ()
instance NFData Bounded.Verdict where rnf !_ = ()
+-- partial helper for benches: configs here are hardcoded valid, so a
+-- 'Left' would be a bench-suite bug.
+ok :: Either e a -> a
+ok (Right x) = x
+ok (Left _) = error "bench: invalid config"
+
main :: IO ()
main = defaultMain [
update
, decide
, stream
, twosample
+ , bern_update
+ , bern_stream
]
update :: Benchmark
update =
- let !cfg_f = Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
- !cfg_a = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive
- !cfg_o = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ let !cfg_f = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
+ !cfg_a = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive)
+ !cfg_o = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
!st_f = Bounded.initial cfg_f
!st_a = Bounded.initial cfg_a
!st_o = Bounded.initial cfg_o
!x = 0.7
in bgroup "Bounded.update (one step)" [
- bench "fixed" $ nf (Bounded.update cfg_f st_f) x
+ bench "fixed" $ nf (Bounded.update cfg_f st_f) x
, bench "adaptive" $ nf (Bounded.update cfg_a st_a) x
- , bench "newton" $ nf (Bounded.update cfg_o st_o) x
+ , bench "newton" $ nf (Bounded.update cfg_o st_o) x
]
decide :: Benchmark
decide =
- let !cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ let !cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
!st = Bounded.initial cfg
in bgroup "Bounded.decide" [
bench "initial state" $ nf (Bounded.decide cfg) st
@@ -49,25 +59,52 @@ decide =
stream :: Benchmark
stream =
let !xs = force (take 1000 (cycle [0.3, 0.7]))
- !cfg_f = Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
- !cfg_a = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive
- !cfg_o = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ !cfg_f = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
+ !cfg_a = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive)
+ !cfg_o = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
run_m cfg = foldl' (Bounded.update cfg) (Bounded.initial cfg)
in bgroup "Bounded.update (1000-sample fold)" [
- bench "fixed" $ nf (run_m cfg_f) xs
+ bench "fixed" $ nf (run_m cfg_f) xs
, bench "adaptive" $ nf (run_m cfg_a) xs
- , bench "newton" $ nf (run_m cfg_o) xs
+ , bench "newton" $ nf (run_m cfg_o) xs
]
twosample :: Benchmark
twosample =
let !ps = force (take 1000 (cycle [(0.3, 0.7), (0.7, 0.3)]))
- !cfg_f = P.config 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
- !cfg_a = P.config 0.0 1.0 1.0e-3 Bounded.Adaptive
- !cfg_o = P.config 0.0 1.0 1.0e-3 Bounded.Newton
+ !cfg_f = ok (P.config 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
+ !cfg_a = ok (P.config 0.0 1.0 1.0e-3 Bounded.Adaptive)
+ !cfg_o = ok (P.config 0.0 1.0 1.0e-3 Bounded.Newton)
run_t cfg = foldl' (P.update cfg) (P.initial cfg)
in bgroup "Paired.update (1000-sample fold)" [
- bench "fixed" $ nf (run_t cfg_f) ps
+ bench "fixed" $ nf (run_t cfg_f) ps
, bench "adaptive" $ nf (run_t cfg_a) ps
- , bench "newton" $ nf (run_t cfg_o) ps
+ , bench "newton" $ nf (run_t cfg_o) ps
+ ]
+
+bern_update :: Benchmark
+bern_update =
+ let !cfg_f = ok (Bern.config 0.05 1.0e-3 (Bern.Fixed 5.0))
+ !cfg_a = ok (Bern.config 0.05 1.0e-3 Bern.Adaptive)
+ !cfg_o = ok (Bern.config 0.05 1.0e-3 Bern.Newton)
+ !st_f = Bern.initial cfg_f
+ !st_a = Bern.initial cfg_a
+ !st_o = Bern.initial cfg_o
+ in bgroup "Bernoulli.update (one step)" [
+ bench "fixed" $ nf (Bern.update cfg_f st_f) True
+ , bench "adaptive" $ nf (Bern.update cfg_a st_a) True
+ , bench "newton" $ nf (Bern.update cfg_o st_o) True
+ ]
+
+bern_stream :: Benchmark
+bern_stream =
+ let !xs = force (take 1000 (cycle [True, False]))
+ !cfg_f = ok (Bern.config 0.05 1.0e-3 (Bern.Fixed 5.0))
+ !cfg_a = ok (Bern.config 0.05 1.0e-3 Bern.Adaptive)
+ !cfg_o = ok (Bern.config 0.05 1.0e-3 Bern.Newton)
+ run_b cfg = foldl' (Bern.update cfg) (Bern.initial cfg)
+ in bgroup "Bernoulli.update (1000-sample fold)" [
+ bench "fixed" $ nf (run_b cfg_f) xs
+ , bench "adaptive" $ nf (run_b cfg_a) xs
+ , bench "newton" $ nf (run_b cfg_o) xs
]
diff --git a/bench/Weight.hs b/bench/Weight.hs
@@ -4,14 +4,21 @@
module Main where
import Control.DeepSeq
+import qualified Numeric.Eproc.Bernoulli as Bern
import qualified Numeric.Eproc.Bounded as Bounded
import qualified Numeric.Eproc.Paired as P
import Weigh
instance NFData Bounded.State where rnf !_ = ()
-instance NFData P.State where rnf !_ = ()
+instance NFData P.State where rnf !_ = ()
+instance NFData Bern.State where rnf !_ = ()
instance NFData Bounded.Verdict where rnf !_ = ()
+-- partial helper for benches: configs here are hardcoded valid.
+ok :: Either e a -> a
+ok (Right x) = x
+ok (Left _) = error "weigh: invalid config"
+
-- note that 'weigh' doesn't work properly in a repl
main :: IO ()
main = mainWith $ do
@@ -19,23 +26,25 @@ main = mainWith $ do
decide
stream
twosample
+ bern_update
+ bern_stream
update :: Weigh ()
update =
- let !cfg_f = Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
- !cfg_a = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive
- !cfg_o = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ let !cfg_f = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
+ !cfg_a = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive)
+ !cfg_o = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
!st_f = Bounded.initial cfg_f
!st_a = Bounded.initial cfg_a
!st_o = Bounded.initial cfg_o
in wgroup "Bounded.update (one step)" $ do
- func "fixed" (Bounded.update cfg_f st_f) 0.7
+ func "fixed" (Bounded.update cfg_f st_f) 0.7
func "adaptive" (Bounded.update cfg_a st_a) 0.7
- func "newton" (Bounded.update cfg_o st_o) 0.7
+ func "newton" (Bounded.update cfg_o st_o) 0.7
decide :: Weigh ()
decide =
- let !cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ let !cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
!st = Bounded.initial cfg
in wgroup "Bounded.decide" $ do
func "initial state" (Bounded.decide cfg) st
@@ -43,23 +52,48 @@ decide =
stream :: Weigh ()
stream =
let !xs = force (take 1000 (cycle [0.3, 0.7]))
- !cfg_f = Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
- !cfg_a = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive
- !cfg_o = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ !cfg_f = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
+ !cfg_a = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive)
+ !cfg_o = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
run_m cfg = foldl' (Bounded.update cfg) (Bounded.initial cfg)
in wgroup "Bounded.update (1000-sample fold)" $ do
- func "fixed" (run_m cfg_f) xs
+ func "fixed" (run_m cfg_f) xs
func "adaptive" (run_m cfg_a) xs
- func "newton" (run_m cfg_o) xs
+ func "newton" (run_m cfg_o) xs
twosample :: Weigh ()
twosample =
let !ps = force (take 1000 (cycle [(0.3, 0.7), (0.7, 0.3)]))
- !cfg_f = P.config 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
- !cfg_a = P.config 0.0 1.0 1.0e-3 Bounded.Adaptive
- !cfg_o = P.config 0.0 1.0 1.0e-3 Bounded.Newton
+ !cfg_f = ok (P.config 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
+ !cfg_a = ok (P.config 0.0 1.0 1.0e-3 Bounded.Adaptive)
+ !cfg_o = ok (P.config 0.0 1.0 1.0e-3 Bounded.Newton)
run_t cfg = foldl' (P.update cfg) (P.initial cfg)
in wgroup "Paired.update (1000-sample fold)" $ do
- func "fixed" (run_t cfg_f) ps
+ func "fixed" (run_t cfg_f) ps
func "adaptive" (run_t cfg_a) ps
- func "newton" (run_t cfg_o) ps
+ func "newton" (run_t cfg_o) ps
+
+bern_update :: Weigh ()
+bern_update =
+ let !cfg_f = ok (Bern.config 0.05 1.0e-3 (Bern.Fixed 5.0))
+ !cfg_a = ok (Bern.config 0.05 1.0e-3 Bern.Adaptive)
+ !cfg_o = ok (Bern.config 0.05 1.0e-3 Bern.Newton)
+ !st_f = Bern.initial cfg_f
+ !st_a = Bern.initial cfg_a
+ !st_o = Bern.initial cfg_o
+ in wgroup "Bernoulli.update (one step)" $ do
+ func "fixed" (Bern.update cfg_f st_f) True
+ func "adaptive" (Bern.update cfg_a st_a) True
+ func "newton" (Bern.update cfg_o st_o) True
+
+bern_stream :: Weigh ()
+bern_stream =
+ let !xs = force (take 1000 (cycle [True, False]))
+ !cfg_f = ok (Bern.config 0.05 1.0e-3 (Bern.Fixed 5.0))
+ !cfg_a = ok (Bern.config 0.05 1.0e-3 Bern.Adaptive)
+ !cfg_o = ok (Bern.config 0.05 1.0e-3 Bern.Newton)
+ run_b cfg = foldl' (Bern.update cfg) (Bern.initial cfg)
+ in wgroup "Bernoulli.update (1000-sample fold)" $ do
+ func "fixed" (run_b cfg_f) xs
+ func "adaptive" (run_b cfg_a) xs
+ func "newton" (run_b cfg_o) xs
diff --git a/lib/Numeric/Eproc/Bernoulli.hs b/lib/Numeric/Eproc/Bernoulli.hs
@@ -10,8 +10,15 @@
--
-- One-sided Bernoulli rate anytime-valid test.
--
--- For samples @x_t@ in @{0, 1}@, tests @H_0: E[x] <= p_0@ against
--- @H_1: E[x] > p_0@.
+-- For samples @x_t@ in @{0, 1}@, tests
+--
+-- @H_0: E[x_t | F_{t-1}] <= p_0 for all t@
+--
+-- against @H_1: E[x_t | F_{t-1}] > p_0@ (at some @t@). Here
+-- @F_{t-1}@ is the filtration generated by everything observed
+-- strictly before time @t@; the conditional form is what anytime
+-- validity actually requires. For i.i.d. samples this reduces to
+-- the usual marginal statement @E[x] <= p_0@.
--
-- A single wealth process is run:
--
@@ -24,7 +31,10 @@
-- a nonnegative supermartingale, so by Ville's inequality the
-- probability of @W_n@ ever crossing @1 \/ alpha@ is at most
-- @alpha@, regardless of when the user decides to stop streaming
--- samples.
+-- samples. Rejection is /latched/ in the running state: once the
+-- wealth has crossed threshold, 'decide' continues to return
+-- 'Reject' even if subsequent observations drive the current
+-- wealth back below threshold.
--
-- Unlike "Numeric.Eproc.Bounded", the alternative here is one-sided,
-- so a single wealth process suffices and no Bonferroni adjustment
@@ -35,7 +45,7 @@
-- Test @H_0: E[x] <= 0.05@ at level @alpha = 1e-3@ against a stream
-- with empirical rate @~0.5@:
--
--- >>> let cfg = config 1.0e-3 0.05 Newton
+-- >>> let Right cfg = config 0.05 1.0e-3 Newton
-- >>> let xs = take 200 (cycle [True, False])
-- >>> decide cfg (foldl' (update cfg) (initial cfg) xs)
-- Reject
@@ -45,6 +55,7 @@ module Numeric.Eproc.Bernoulli (
Config
, State
, Verdict(..)
+ , ConfigError(..)
-- * Bettor strategies
, Bettor(..)
@@ -62,7 +73,10 @@ module Numeric.Eproc.Bernoulli (
, samples
) where
-import Numeric.Eproc.Common (Bettor(..), Verdict(..))
+import Numeric.Eproc.Common (
+ Bettor(..), Verdict(..), ConfigError(..)
+ , BetState, init_bet, bet_lambda, step_bet
+ )
-- types ----------------------------------------------------------------------
@@ -70,19 +84,6 @@ import Numeric.Eproc.Common (Bettor(..), Verdict(..))
-- "Numeric.Eproc.Common" is @x_t - p_0@; the safe-bet ceiling
-- @lambda_max@ is derived from @p_0@ (see 'config').
--- bettor state. one constructor per 'Bettor' alternative; the
--- constructor used in a given 'State' matches the 'Bettor' chosen in
--- the enclosing 'Config'.
-data BetState =
- SFixed
- | SAdaptive
- {-# UNPACK #-} !Double -- sum of z (centred observation)
- {-# UNPACK #-} !Double -- sum of z^2 (for online variance)
- {-# UNPACK #-} !Int -- count
- | SNewton
- {-# UNPACK #-} !Double -- current bet lambda
- {-# UNPACK #-} !Double -- running sum of per-step squared gradients
-
-- | Bernoulli rate test configuration. Build with 'config'.
--
-- Carries the bettor strategy, the baseline rate, the significance
@@ -104,69 +105,18 @@ data Config = Config {
-- | Streaming test state. Construct with 'initial' and fold
-- observations through 'update'.
--
--- Carries the sample count, running log-wealth, and whatever
--- per-step state the chosen 'Bettor' needs.
+-- Carries the sample count, current and supremum-so-far running
+-- log-wealth, and whatever per-step state the chosen 'Bettor'
+-- needs. The supremum field is what 'decide' tests against the
+-- rejection threshold; this is the supremum-style event Ville's
+-- inequality actually bounds.
data State = State {
- st_n :: {-# UNPACK #-} !Int -- ^ sample count
- , st_log_w :: {-# UNPACK #-} !Double -- ^ running log-wealth
- , st_bet :: !BetState -- ^ bettor state
+ st_n :: {-# UNPACK #-} !Int -- ^ sample count
+ , st_log_w :: {-# UNPACK #-} !Double -- ^ running log-wealth
+ , st_max_log_w :: {-# UNPACK #-} !Double -- ^ sup log-wealth so far
+ , st_bet :: !BetState -- ^ bettor state
}
--- internal -------------------------------------------------------------------
-
--- per-bettor initial state.
-init_bet :: Bettor -> BetState
-init_bet b = case b of
- Fixed _ -> SFixed
- Adaptive -> SAdaptive 0 0 0
- Newton -> SNewton 0 1.0e-6 -- small acc seed avoids div-by-zero
-{-# INLINE init_bet #-}
-
--- compute the next bet 'lambda' from the bettor and its current
--- state. for Adaptive we form a Kelly-style plug-in from the running
--- sample mean and variance; for Newton the bet is just the last
--- lambda chosen by the Newton step (updated during 'step_bet').
-bet_lambda :: Bettor -> Double -> BetState -> Double
-bet_lambda b !lam_max !s = case b of
- Fixed lam -> lam
- Adaptive -> case s of
- SAdaptive !sm !sm2 !n
- | n == 0 -> 0
- | otherwise ->
- let !nd = fromIntegral n
- !mu = sm / nd
- !mu2 = mu * mu
- !var = max 0 (sm2 / nd - mu2)
- !den = var + mu2
- !raw = if den == 0 then 0 else mu / den
- in max 0 (min lam_max raw)
- _ -> 0
- Newton -> case s of
- SNewton !lam _ -> lam
- _ -> 0
-{-# INLINE bet_lambda #-}
-
--- update bettor state with newly observed centred value 'z'. for
--- Adaptive this is just accumulating sums; for Newton we take one
--- Newton step on the per-step log-wealth loss '-log(1 + lambda * z)',
--- accumulating squared gradients for adaptive scaling.
-step_bet :: Bettor -> Double -> BetState -> Double -> BetState
-step_bet b !lam_max !s !z = case b of
- Fixed _ -> SFixed
- Adaptive -> case s of
- SAdaptive !sm !sm2 !n -> SAdaptive (sm + z) (sm2 + z * z) (n + 1)
- _ -> SAdaptive z (z * z) 1
- Newton -> case s of
- SNewton !lam !acc ->
- let !denom = 1 + lam * z
- !g = if denom == 0 then 0 else negate z / denom
- !acc' = acc + g * g
- !lam' = lam - g / acc'
- !clp = max 0 (min lam_max lam')
- in SNewton clp acc'
- _ -> SNewton 0 1.0e-6
-{-# INLINE step_bet #-}
-
-- construction ---------------------------------------------------------------
-- | Build a 'Config' for the Bernoulli rate test.
@@ -177,39 +127,42 @@ step_bet b !lam_max !s !z = case b of
-- requires @lambda <= 1 \/ p_0@; the ceiling stored is half this
-- to leave numerical margin -- the WSR safety recommendation.
--
--- @p_0@ must lie strictly in @(0, 1)@ and @alpha@ strictly in
--- @(0, 1)@. The degenerate case @p_0 = 0@ would make @lambda_max@
--- infinite (any divergence would reject immediately and the test
--- becomes uninteresting); the caller is expected to pass a small
--- positive baseline.
+-- Returns 'Left' with a 'ConfigError' on inputs that would leave
+-- the mathematical regime: @p_0@ outside @(0, 1)@ (the degenerate
+-- case @p_0 = 0@ would make @lambda_max@ infinite, and @p_0 = 1@
+-- leaves no room for an alternative), or @alpha@ outside @(0, 1)@.
--
--- >>> let cfg = config 1.0e-3 0.05 Newton
+-- >>> let Right cfg = config 0.05 1.0e-3 Newton
config
- :: Double -- ^ significance level @alpha@, in @(0, 1)@
- -> Double -- ^ baseline rate @p_0@, in @(0, 1)@
+ :: Double -- ^ baseline rate @p_0@, in @(0, 1)@
+ -> Double -- ^ significance level @alpha@, in @(0, 1)@
-> Bettor -- ^ bettor strategy
- -> Config
-config !alpha !p0 !b = Config {
- cfg_bettor = b
- , cfg_lam_max = 0.5 / p0
- , cfg_p0 = p0
- , cfg_alpha = alpha
- , cfg_log_thresh = log (1 / alpha)
- }
+ -> Either ConfigError Config
+config !p0 !alpha !b
+ | not (p0 > 0 && p0 < 1) = Left (InvalidBaselineRate p0)
+ | not (alpha > 0 && alpha < 1) = Left (InvalidAlpha alpha)
+ | otherwise = Right Config {
+ cfg_bettor = b
+ , cfg_lam_max = 0.5 / p0
+ , cfg_p0 = p0
+ , cfg_alpha = alpha
+ , cfg_log_thresh = log (1 / alpha)
+ }
{-# INLINE config #-}
-- | The initial 'State' for a fresh streaming test.
--
--- Log-wealth starts at @0@ (i.e., wealth @1@) and the bettor
--- starts in the per-strategy initial state appropriate for the
--- 'Bettor' chosen in the 'Config'.
+-- Both log-wealth fields start at @0@ (i.e., wealth @1@) and the
+-- bettor starts in the per-strategy initial state appropriate
+-- for the 'Bettor' chosen in the 'Config'.
--
-- >>> let s0 = initial cfg
initial :: Config -> State
initial Config{..} = State {
- st_n = 0
- , st_log_w = 0
- , st_bet = init_bet cfg_bettor
+ st_n = 0
+ , st_log_w = 0
+ , st_max_log_w = 0
+ , st_bet = init_bet cfg_bettor
}
{-# INLINE initial #-}
@@ -227,7 +180,12 @@ initial Config{..} = State {
--
-- @log_w' = log_w + log (1 + lambda * z)@
--
--- and then steps the bettor state given the newly observed @z@.
+-- updates the running supremum log-wealth, then steps the bettor
+-- state given the newly observed @z@.
+--
+-- /Precondition/: @True@ and @False@ both /must/ be admissible
+-- under the test (this holds vacuously for the @{0, 1}@ support).
+-- The function is total.
--
-- >>> let s1 = update cfg s0 True
update :: Config -> State -> Bool -> State
@@ -237,41 +195,42 @@ update Config{..} State{..} !x =
!lam = bet_lambda cfg_bettor cfg_lam_max st_bet
!fac = 1 + lam * z
!logw' = st_log_w + log fac
+ !maxw' = max st_max_log_w logw'
!s' = step_bet cfg_bettor cfg_lam_max st_bet z
- in State (st_n + 1) logw' s'
+ in State (st_n + 1) logw' maxw' s'
{-# INLINE update #-}
-- | Compute the current 'Verdict' from the running 'State'.
--
--- 'Reject' iff log-wealth has crossed the threshold
+-- 'Reject' iff log-wealth has /ever/ crossed the threshold
-- @log(1 \/ alpha)@; equivalently, wealth has exceeded
--- @1 \/ alpha@. Under @H_0@, by Ville's inequality, the
--- probability of this ever happening is at most @alpha@ -- and
--- crucially this bound holds at /every/ sample size
--- simultaneously, so the user is free to peek at the verdict as
--- often as they like and stop on the first 'Reject'.
+-- @1 \/ alpha@ at some point in the stream so far. Under @H_0@,
+-- by Ville's inequality, the probability of this ever happening
+-- is at most @alpha@ -- and crucially this bound holds at /every/
+-- sample size simultaneously, so the user is free to peek at the
+-- verdict as often as they like and stop on the first 'Reject'.
--
-- >>> decide cfg s0
-- Continue
decide :: Config -> State -> Verdict
decide Config{..} State{..}
- | st_log_w >= cfg_log_thresh = Reject
- | otherwise = Continue
+ | st_max_log_w >= cfg_log_thresh = Reject
+ | otherwise = Continue
{-# INLINE decide #-}
-- inspection -----------------------------------------------------------------
--- | The current log-wealth.
+-- | The supremum-so-far log-wealth, across all sample counts up to
+-- the current one.
--
--- This is the natural \"test statistic\": it is monotone (in
--- expectation under @H_1@) in the evidence against @H_0@
--- accumulated so far, and the test rejects exactly when it crosses
--- @log(1 \/ alpha)@.
+-- This is the natural \"test statistic\": it is monotone
+-- nondecreasing in the sample count, and 'decide' rejects exactly
+-- when it crosses @log(1 \/ alpha)@.
--
-- >>> log_wealth s0
-- 0.0
log_wealth :: State -> Double
-log_wealth = st_log_w
+log_wealth = st_max_log_w
{-# INLINE log_wealth #-}
-- | The number of samples consumed so far.
diff --git a/lib/Numeric/Eproc/Bounded.hs b/lib/Numeric/Eproc/Bounded.hs
@@ -1,6 +1,5 @@
{-# OPTIONS_HADDOCK prune #-}
{-# LANGUAGE BangPatterns #-}
-{-# LANGUAGE MagicHash #-}
{-# LANGUAGE RecordWildCards #-}
-- |
@@ -11,28 +10,41 @@
--
-- Two-sided bounded-mean anytime-valid test.
--
--- For samples @x_t@ in @[lo, hi]@, tests @H_0: E[x] = m@ against
--- @H_1: E[x] /= m@.
+-- For samples @x_t@ in @[lo, hi]@, tests
+--
+-- @H_0: E[x_t | F_{t-1}] = m for all t@
+--
+-- against the negation. Here @F_{t-1}@ is the filtration generated
+-- by everything observed strictly before time @t@; the conditional
+-- form is what anytime validity actually requires. For i.i.d.
+-- samples this reduces to the usual marginal statement
+-- @E[x] = m@; for adaptively-collected or otherwise non-i.i.d.
+-- streams the conditional statement is the right thing to think
+-- about.
--
-- Internally two one-sided e-processes are run in parallel: a
-- /positive-direction/ process betting against the alternative
--- @E[x] > m@ (using centred observations @z = x - m@), and a
--- /negative-direction/ process betting against @E[x] < m@ (using
--- @-z@). Each maintains its own log-wealth and bettor state. The
--- test rejects when either side's wealth crosses @2 \/ alpha@; the
--- factor of 2 is the Bonferroni adjustment for the two-sided union.
+-- @E[x_t | F_{t-1}] > m@ (using centred observations @z = x - m@),
+-- and a /negative-direction/ process betting against
+-- @E[x_t | F_{t-1}] < m@ (using @-z@). Each maintains its own
+-- log-wealth and bettor state. The test rejects when /either/
+-- side's wealth has /ever/ crossed @2 \/ alpha@; the factor of 2
+-- is the Bonferroni adjustment for the two-sided union.
--
-- The test is /anytime-valid/: under @H_0@ the wealth process is a
-- nonnegative supermartingale, so by Ville's inequality the
--- probability of ever crossing the threshold is at most @alpha@,
+-- probability of /ever/ crossing the threshold is at most @alpha@,
-- regardless of when the user decides to stop streaming samples.
+-- Rejection is /latched/ in the running state -- once a side has
+-- crossed threshold, 'decide' continues to return 'Reject' even if
+-- the current log-wealth has since dropped back below threshold.
--
-- == Example
--
-- Test @H_0: E[x] = 0.5@ for @x@ in @[0, 1]@ at level @alpha = 1e-3@
-- against a stream with empirical mean @0.8@:
--
--- >>> let cfg = config 0.5 0.0 1.0 1.0e-3 Newton
+-- >>> let Right cfg = config 0.5 0.0 1.0 1.0e-3 Newton
-- >>> let xs = concat (replicate 30 [1, 1, 0, 1, 1, 0, 1, 1, 1, 1])
-- >>> decide cfg (foldl' (update cfg) (initial cfg) xs)
-- Reject
@@ -42,6 +54,7 @@ module Numeric.Eproc.Bounded (
Config
, State
, Verdict(..)
+ , ConfigError(..)
-- * Bettor strategies
, Bettor(..)
@@ -59,8 +72,10 @@ module Numeric.Eproc.Bounded (
, samples
) where
-import GHC.Exts (Double(D#))
-import Numeric.Eproc.Common (Bettor(..), Verdict(..))
+import Numeric.Eproc.Common (
+ Bettor(..), Verdict(..), ConfigError(..)
+ , BetState, init_bet, bet_lambda, step_bet
+ )
-- types ----------------------------------------------------------------------
@@ -69,19 +84,6 @@ import Numeric.Eproc.Common (Bettor(..), Verdict(..))
-- ceilings @lambda_max@ are derived from the sample bounds (see
-- 'config').
--- per-direction bettor state. one constructor per 'Bettor' alternative;
--- the constructor used in a given 'State' matches the 'Bettor' chosen
--- in the enclosing 'Config'.
-data BetState =
- SFixed
- | SAdaptive
- {-# UNPACK #-} !Double -- sum of z (centred observation)
- {-# UNPACK #-} !Double -- sum of z^2 (for online variance)
- {-# UNPACK #-} !Int -- count
- | SNewton
- {-# UNPACK #-} !Double -- current bet lambda
- {-# UNPACK #-} !Double -- running sum of per-step squared gradients
-
-- | Bounded-mean test configuration. Build with 'config'.
--
-- Carries the bettor strategy, the null mean, the significance
@@ -107,85 +109,22 @@ data Config = Config {
-- observations through 'update'.
--
-- The two log-wealth fields track the running log-wealth of the
--- positive- and negative-direction e-processes separately;
--- 'decide' compares each to the threshold and 'log_wealth' returns
--- the larger of the two. The per-direction bettor states carry
--- whatever the chosen 'Bettor' needs (running sums, current bet,
--- etc.).
+-- positive- and negative-direction e-processes separately; the
+-- two /maximum/ log-wealth fields latch the supremum so far on
+-- each side, so 'decide' tests the supremum-style event Ville's
+-- inequality actually bounds. The per-direction bettor states
+-- carry whatever the chosen 'Bettor' needs (running sums, current
+-- bet, etc.).
data State = State {
- st_n :: {-# UNPACK #-} !Int -- ^ sample count
- , st_log_w_pos :: {-# UNPACK #-} !Double -- ^ log-wealth, pos-dir process
- , st_log_w_neg :: {-# UNPACK #-} !Double -- ^ log-wealth, neg-dir process
- , st_bet_pos :: !BetState -- ^ bettor state, pos-direction
- , st_bet_neg :: !BetState -- ^ bettor state, neg-direction
+ st_n :: {-# UNPACK #-} !Int -- ^ sample count
+ , st_log_w_pos :: {-# UNPACK #-} !Double -- ^ log-wealth, pos
+ , st_log_w_neg :: {-# UNPACK #-} !Double -- ^ log-wealth, neg
+ , st_max_log_w_pos :: {-# UNPACK #-} !Double -- ^ sup log-wealth, pos
+ , st_max_log_w_neg :: {-# UNPACK #-} !Double -- ^ sup log-wealth, neg
+ , st_bet_pos :: !BetState -- ^ bettor state, pos
+ , st_bet_neg :: !BetState -- ^ bettor state, neg
}
--- internal -------------------------------------------------------------------
-
--- floor for the wealth factor before taking a log; keeps the running
--- log-wealth finite when a step pushes the factor to (or below) zero.
--- NB. written via MagicHash because the fractional literal '1.0e-300'
--- compiles as 'fromRational (1.0e-300 :: Rational)', and GHC does
--- not constant-fold the conversion -- leaving a per-step
--- '$wrationalToDouble' call in the worker.
-tiny :: Double
-tiny = D# 1.0e-300##
-{-# INLINE tiny #-}
-
--- per-bettor initial state.
-init_bet :: Bettor -> BetState
-init_bet b = case b of
- Fixed _ -> SFixed
- Adaptive -> SAdaptive 0 0 0
- Newton -> SNewton 0 1.0e-6 -- small acc seed avoids div-by-zero
-{-# INLINE init_bet #-}
-
--- compute the next bet 'lambda' from the bettor and its current
--- state; 'lam_max' is the direction-specific safety bound. for
--- Adaptive we form a Kelly-style plug-in from the running sample
--- mean and variance; for Newton the bet is just the last lambda
--- chosen by the Newton step (updated during 'step_bet').
-bet_lambda :: Bettor -> Double -> BetState -> Double
-bet_lambda b !lam_max !s = case b of
- Fixed lam -> lam
- Adaptive -> case s of
- SAdaptive !sm !sm2 !n
- | n == 0 -> 0
- | otherwise ->
- let !nd = fromIntegral n
- !mu = sm / nd
- !mu2 = mu * mu
- !var = max 0 (sm2 / nd - mu2)
- !den = var + mu2
- !raw = if den == 0 then 0 else mu / den
- in max 0 (min lam_max raw)
- _ -> 0
- Newton -> case s of
- SNewton !lam _ -> lam
- _ -> 0
-{-# INLINE bet_lambda #-}
-
--- update bettor state with newly observed centred value 'z'. for
--- Adaptive this is just accumulating sums; for Newton we take one
--- Newton step on the per-step log-wealth loss '-log(1 + lambda * z)',
--- accumulating squared gradients for adaptive scaling.
-step_bet :: Bettor -> Double -> BetState -> Double -> BetState
-step_bet b !lam_max !s !z = case b of
- Fixed _ -> SFixed
- Adaptive -> case s of
- SAdaptive !sm !sm2 !n -> SAdaptive (sm + z) (sm2 + z * z) (n + 1)
- _ -> SAdaptive z (z * z) 1
- Newton -> case s of
- SNewton !lam !acc ->
- let !denom = 1 + lam * z
- !g = if denom == 0 then 0 else negate z / denom
- !acc' = acc + g * g
- !lam' = lam - g / acc'
- !clp = max 0 (min lam_max lam')
- in SNewton clp acc'
- _ -> SNewton 0 1.0e-6
-{-# INLINE step_bet #-}
-
-- construction ---------------------------------------------------------------
-- | Build a 'Config' for the bounded-mean test.
@@ -209,27 +148,36 @@ step_bet b !lam_max !s !z = case b of
-- @log(2 \/ alpha)@; the 2 is the Bonferroni union-bound
-- adjustment for the two one-sided e-processes.
--
--- >>> let cfg = config 0.5 0.0 1.0 1.0e-3 Newton
+-- Returns 'Left' with a 'ConfigError' on inputs that would leave
+-- the mathematical regime: @alpha@ outside @(0, 1)@, @lo >= hi@,
+-- or @m@ outside the open interval @(lo, hi)@ (strict, to avoid
+-- the safe-bet ceilings dividing by zero).
+--
+-- >>> let Right cfg = config 0.5 0.0 1.0 1.0e-3 Newton
config
:: Double -- ^ null mean @m@
-> Double -- ^ sample lower bound @lo@
-> Double -- ^ sample upper bound @hi@
-> Double -- ^ significance level @alpha@
-> Bettor -- ^ bettor strategy
- -> Config
-config !m !lo !hi !alpha !b = Config {
- cfg_bettor = b
- , cfg_lam_max_pos = 0.5 / (m - lo)
- , cfg_lam_max_neg = 0.5 / (hi - m)
- , cfg_null_mean = m
- , cfg_alpha = alpha
- , cfg_log_thresh = log (2 / alpha)
- }
+ -> Either ConfigError Config
+config !m !lo !hi !alpha !b
+ | not (alpha > 0 && alpha < 1) = Left (InvalidAlpha alpha)
+ | not (lo < hi) = Left (InvalidBounds lo hi)
+ | not (lo < m && m < hi) = Left (InvalidNullMean m lo hi)
+ | otherwise = Right Config {
+ cfg_bettor = b
+ , cfg_lam_max_pos = 0.5 / (m - lo)
+ , cfg_lam_max_neg = 0.5 / (hi - m)
+ , cfg_null_mean = m
+ , cfg_alpha = alpha
+ , cfg_log_thresh = log (2 / alpha)
+ }
{-# INLINE config #-}
-- | The initial 'State' for a fresh streaming test.
--
--- Both directional log-wealths start at @0@ (i.e., wealth @1@) and
+-- All four log-wealth fields start at @0@ (i.e., wealth @1@), and
-- both bettors start in the per-strategy initial state appropriate
-- for the 'Bettor' chosen in the 'Config'.
--
@@ -238,11 +186,13 @@ initial :: Config -> State
initial Config{..} =
let !s0 = init_bet cfg_bettor
in State {
- st_n = 0
- , st_log_w_pos = 0
- , st_log_w_neg = 0
- , st_bet_pos = s0
- , st_bet_neg = s0
+ st_n = 0
+ , st_log_w_pos = 0
+ , st_log_w_neg = 0
+ , st_max_log_w_pos = 0
+ , st_max_log_w_neg = 0
+ , st_bet_pos = s0
+ , st_bet_neg = s0
}
{-# INLINE initial #-}
@@ -256,11 +206,15 @@ initial Config{..} =
--
-- @log_w' = log_w + log (1 + lambda * z)@
--
--- (with the symmetric @-lambda@ for the negative direction), and
--- then steps the bettor states given the newly observed @z@. The
--- per-step wealth factor is floored at a tiny positive value to
--- keep the log finite when a marginal bet drives the factor to (or
--- below) zero.
+-- (with the symmetric @-lambda@ for the negative direction), then
+-- updates the running supremum of log-wealth on each side and
+-- steps the bettor states given the newly observed @z@.
+--
+-- /Precondition/: @x@ must lie in the @[lo, hi]@ interval given
+-- to 'config'. The type-I error guarantee of the test depends on
+-- this. Out-of-range observations can drive the wealth factor
+-- negative, taking the construction out of the supermartingale
+-- regime entirely; the function does not check for this.
--
-- >>> let s1 = update cfg s0 0.7
update :: Config -> State -> Double -> State
@@ -270,46 +224,49 @@ update Config{..} State{..} !x =
!lam_n = bet_lambda cfg_bettor cfg_lam_max_neg st_bet_neg
!fac_p = 1 + lam_p * z
!fac_n = 1 - lam_n * z
- !logw_p = st_log_w_pos + log (max tiny fac_p)
- !logw_n = st_log_w_neg + log (max tiny fac_n)
+ !logw_p = st_log_w_pos + log fac_p
+ !logw_n = st_log_w_neg + log fac_n
+ !maxp = max st_max_log_w_pos logw_p
+ !maxn = max st_max_log_w_neg logw_n
!sp = step_bet cfg_bettor cfg_lam_max_pos st_bet_pos z
!sn = step_bet cfg_bettor cfg_lam_max_neg st_bet_neg (negate z)
- in State (st_n + 1) logw_p logw_n sp sn
+ in State (st_n + 1) logw_p logw_n maxp maxn sp sn
{-# INLINE update #-}
-- | Compute the current 'Verdict' from the running 'State'.
--
--- 'Reject' iff either directional log-wealth has crossed the
--- Bonferroni-adjusted threshold @log(2 \/ alpha)@; equivalently,
--- the wealth process on either side has exceeded @2 \/ alpha@.
--- Under @H_0@, by Ville's inequality, the probability of this ever
--- happening is at most @alpha@ -- and crucially this bound holds
--- at /every/ sample size simultaneously, so the user is free to
--- peek at the verdict as often as they like and stop on the first
--- 'Reject'.
+-- 'Reject' iff either directional log-wealth has /ever/ crossed
+-- the Bonferroni-adjusted threshold @log(2 \/ alpha)@;
+-- equivalently, the wealth process on either side has exceeded
+-- @2 \/ alpha@ at some point in the stream so far. Under @H_0@,
+-- by Ville's inequality, the probability of this ever happening
+-- is at most @alpha@ -- and crucially this bound holds at /every/
+-- sample size simultaneously, so the user is free to peek at the
+-- verdict as often as they like and stop on the first 'Reject'.
--
-- >>> decide cfg s0
-- Continue
decide :: Config -> State -> Verdict
decide Config{..} State{..}
- | st_log_w_pos >= cfg_log_thresh = Reject
- | st_log_w_neg >= cfg_log_thresh = Reject
- | otherwise = Continue
+ | st_max_log_w_pos >= cfg_log_thresh = Reject
+ | st_max_log_w_neg >= cfg_log_thresh = Reject
+ | otherwise = Continue
{-# INLINE decide #-}
-- inspection -----------------------------------------------------------------
--- | The current log-wealth, taken as the maximum of the two
--- directional processes.
+-- | The supremum-so-far log-wealth, taken as the maximum across the
+-- two directional processes and across all sample counts up to
+-- the current one.
--
--- This is the natural \"test statistic\": it is monotone in the
--- evidence against @H_0@ accumulated so far, and the test rejects
--- exactly when it crosses @log(2 \/ alpha)@.
+-- This is the natural \"test statistic\": it is monotone
+-- nondecreasing in the sample count, and 'decide' rejects exactly
+-- when it crosses @log(2 \/ alpha)@.
--
-- >>> log_wealth s0
-- 0.0
log_wealth :: State -> Double
-log_wealth State{..} = max st_log_w_pos st_log_w_neg
+log_wealth State{..} = max st_max_log_w_pos st_max_log_w_neg
{-# INLINE log_wealth #-}
-- | The number of samples consumed so far.
diff --git a/lib/Numeric/Eproc/Common.hs b/lib/Numeric/Eproc/Common.hs
@@ -1,4 +1,5 @@
{-# OPTIONS_HADDOCK prune #-}
+{-# LANGUAGE BangPatterns #-}
-- |
-- Module: Numeric.Eproc.Common
@@ -7,14 +8,27 @@
-- Maintainer: Jared Tobin <jared@ppad.tech>
--
-- Shared vocabulary for the eproc tests: the predictable bettor
--- strategies and the test verdict type. Re-exported from each test
--- module ("Numeric.Eproc.Bounded", "Numeric.Eproc.Paired",
+-- strategies, the test verdict type, and the configuration-error
+-- type. Re-exported from each test module
+-- ("Numeric.Eproc.Bounded", "Numeric.Eproc.Paired",
-- "Numeric.Eproc.Bernoulli"); import this module directly only if
-- you need the types without picking a particular test.
+--
+-- The 'BetState' type and its helpers are internal to the library:
+-- they are exposed here so that 'Numeric.Eproc.Bounded' and
+-- 'Numeric.Eproc.Bernoulli' can share one implementation, not for
+-- direct use.
module Numeric.Eproc.Common (
Bettor(..)
, Verdict(..)
+ , ConfigError(..)
+
+ -- * Internal: shared bettor state
+ , BetState(..)
+ , init_bet
+ , bet_lambda
+ , step_bet
) where
-- | A predictable bettor.
@@ -27,14 +41,14 @@ module Numeric.Eproc.Common (
-- what makes the resulting wealth process a nonnegative
-- supermartingale under @H_0@.
--
--- For 'Adaptive' and 'Newton', a safe-bet ceiling @lambda_max@
--- derived from the test's admissible-observation range is enforced
--- by clipping @lambda@ to @[0, lambda_max]@, so the wealth factor
--- stays nonnegative.
+-- All three bettors enforce a safe-bet ceiling @lambda_max@
+-- derived from the test's admissible-observation range by clipping
+-- @lambda@ to @[0, lambda_max]@; this keeps the per-step wealth
+-- factor nonnegative.
--
--- * 'Fixed' always bets the supplied constant @lambda@. The wager
--- does not respond to observed data; this strategy is useful
--- only as a baseline.
+-- * 'Fixed' bets the supplied constant @lambda@ (clipped to
+-- @[0, lambda_max]@). The wager does not respond to observed
+-- data; this strategy is useful only as a baseline.
--
-- * 'Adaptive' is the aGRAPA (approximate growth-rate adaptive
-- predictable plug-in) bettor of Waudby-Smith & Ramdas (2024).
@@ -44,13 +58,14 @@ module Numeric.Eproc.Common (
-- @[0, lambda_max]@. Fast to compute and competitive in
-- practice.
--
--- * 'Newton' is the online Newton step (ONS) bettor. The per-step
+-- * 'Newton' is the online Newton step (ONS) bettor of
+-- Waudby-Smith & Ramdas (2024, Algorithm 2). The per-step
-- log-wealth loss @-log(1 + lambda * z)@ is convex in @lambda@;
-- ONS performs one Newton step per observation, accumulating
--- squared gradients to scale the update. Achieves logarithmic
--- regret against the best constant bet in hindsight and is in
--- practice the strongest of the three bettors under most signal
--- regimes.
+-- squared gradients to scale the update by a fixed learning
+-- rate @2 \/ (2 - log 3)@. Achieves logarithmic regret against
+-- the best constant bet in hindsight and is in practice the
+-- strongest of the three bettors under most signal regimes.
data Bettor =
Fixed {-# UNPACK #-} !Double
| Adaptive
@@ -59,12 +74,111 @@ data Bettor =
-- | Test outcome at the current sample count.
--
--- 'Reject' means the wealth process has crossed the rejection
--- threshold, so @H_0@ is rejected at level @alpha@. 'Continue'
--- means there is not yet enough evidence; collect more samples
--- (or stop and report no rejection -- the type-I error guarantee
--- holds for /any/ stopping rule).
+-- 'Reject' means the wealth process has /ever/ crossed the
+-- rejection threshold, so @H_0@ is rejected at level @alpha@.
+-- Once a state has rejected it stays rejected, even if subsequent
+-- observations drive the current wealth back below threshold;
+-- this is the supremum-style guarantee that Ville's inequality
+-- actually delivers. 'Continue' means there is not yet enough
+-- evidence; collect more samples (or stop and report no
+-- rejection -- the type-I error guarantee holds for /any/
+-- stopping rule).
data Verdict =
Reject
| Continue
deriving (Eq, Show)
+
+-- | Reasons that a test-configuration smart constructor can reject
+-- its inputs. Returned by 'Numeric.Eproc.Bounded.config',
+-- 'Numeric.Eproc.Bernoulli.config', and
+-- 'Numeric.Eproc.Paired.config'.
+data ConfigError =
+ -- | significance level outside @(0, 1)@
+ InvalidAlpha {-# UNPACK #-} !Double
+ -- | sample bounds violate @lo < hi@
+ | InvalidBounds {-# UNPACK #-} !Double {-# UNPACK #-} !Double
+ -- | null mean outside @(lo, hi)@ (strict, to avoid div-by-zero
+ -- in the safe-bet ceilings)
+ | InvalidNullMean
+ {-# UNPACK #-} !Double -- m
+ {-# UNPACK #-} !Double -- lo
+ {-# UNPACK #-} !Double -- hi
+ -- | baseline rate outside @(0, 1)@
+ | InvalidBaselineRate {-# UNPACK #-} !Double
+ deriving (Eq, Show)
+
+-- | Per-bettor state. One constructor per 'Bettor' alternative; the
+-- constructor used in any given state matches the 'Bettor' chosen
+-- in the enclosing 'Config'.
+--
+-- Internal: exposed only so that the per-test 'State' types in
+-- "Numeric.Eproc.Bounded" and "Numeric.Eproc.Bernoulli" can share
+-- one implementation.
+data BetState =
+ SFixed
+ | SAdaptive
+ {-# UNPACK #-} !Double -- sum of z (centred observation)
+ {-# UNPACK #-} !Double -- sum of z^2 (for online variance)
+ {-# UNPACK #-} !Int -- count
+ | SNewton
+ {-# UNPACK #-} !Double -- current bet lambda
+ {-# UNPACK #-} !Double -- running sum of per-step squared gradients
+
+-- | Per-bettor initial state.
+init_bet :: Bettor -> BetState
+init_bet b = case b of
+ Fixed _ -> SFixed
+ Adaptive -> SAdaptive 0 0 0
+ Newton -> SNewton 0 1.0e-6 -- small acc seed avoids div-by-zero
+{-# INLINE init_bet #-}
+
+-- | WSR (2024) Algorithm 2 ONS learning rate, @2 \/ (2 - log 3)@.
+ons_lr :: Double
+ons_lr = 2 / (2 - log 3)
+{-# INLINE ons_lr #-}
+
+-- | Compute the next bet 'lambda' from the bettor and its current
+-- state; 'lam_max' is the direction-specific safety bound. All
+-- strategies clip the result to @[0, lam_max]@ so the wealth
+-- factor stays nonnegative.
+bet_lambda :: Bettor -> Double -> BetState -> Double
+bet_lambda b !lam_max !s = case b of
+ Fixed lam -> max 0 (min lam_max lam)
+ Adaptive -> case s of
+ SAdaptive !sm !sm2 !n
+ | n == 0 -> 0
+ | otherwise ->
+ let !nd = fromIntegral n
+ !mu = sm / nd
+ !mu2 = mu * mu
+ !var = max 0 (sm2 / nd - mu2)
+ !den = var + mu2
+ !raw = if den == 0 then 0 else mu / den
+ in max 0 (min lam_max raw)
+ _ -> 0
+ Newton -> case s of
+ SNewton !lam _ -> lam
+ _ -> 0
+{-# INLINE bet_lambda #-}
+
+-- | Update bettor state with newly observed centred value 'z'. For
+-- 'Adaptive' this is just accumulating sums; for 'Newton' we take
+-- one online Newton step (with the WSR learning rate) on the
+-- per-step log-wealth loss @-log(1 + lambda * z)@, accumulating
+-- squared gradients for adaptive scaling.
+step_bet :: Bettor -> Double -> BetState -> Double -> BetState
+step_bet b !lam_max !s !z = case b of
+ Fixed _ -> SFixed
+ Adaptive -> case s of
+ SAdaptive !sm !sm2 !n -> SAdaptive (sm + z) (sm2 + z * z) (n + 1)
+ _ -> SAdaptive z (z * z) 1
+ Newton -> case s of
+ SNewton !lam !acc ->
+ let !denom = 1 + lam * z
+ !g = if denom == 0 then 0 else negate z / denom
+ !acc' = acc + g * g
+ !lam' = lam - ons_lr * g / acc'
+ !clp = max 0 (min lam_max lam')
+ in SNewton clp acc'
+ _ -> SNewton 0 1.0e-6
+{-# INLINE step_bet #-}
diff --git a/lib/Numeric/Eproc/Paired.hs b/lib/Numeric/Eproc/Paired.hs
@@ -10,14 +10,22 @@
-- Paired two-sample anytime-valid mean-equality test.
--
-- For paired observations @(a_t, b_t)@ where both samples lie in
--- @[lo, hi]@, tests @H_0: E[a] = E[b]@ against
--- @H_1: E[a] /= E[b]@.
+-- @[lo, hi]@, tests
--
--- The reduction is straightforward: under the null, the differences
--- @d_t = a_t - b_t@ have mean zero, and differences of @[lo, hi]@
--- values lie in @[lo - hi, hi - lo]@. So the paired test is just
--- the bounded-mean test ("Numeric.Eproc.Bounded") on @d_t@ with
--- null mean @0@ and sample bounds @[lo - hi, hi - lo]@.
+-- @H_0: E[a_t - b_t | F_{t-1}] = 0 for all t@
+--
+-- against the negation. Here @F_{t-1}@ is the filtration generated
+-- by everything observed strictly before time @t@; the conditional
+-- form is what anytime validity actually requires. For i.i.d. pairs
+-- this reduces to the usual marginal statement @E[a] = E[b]@; for
+-- adaptively-collected or otherwise non-i.i.d. streams the
+-- conditional statement is the right thing to think about.
+--
+-- The reduction is straightforward: under @H_0@, the differences
+-- @d_t = a_t - b_t@ have (conditional) mean zero, and differences
+-- of @[lo, hi]@ values lie in @[lo - hi, hi - lo]@. So the paired
+-- test is just the bounded-mean test ("Numeric.Eproc.Bounded") on
+-- @d_t@ with null mean @0@ and sample bounds @[lo - hi, hi - lo]@.
--
-- Pairing is required: independent two-sample testing without
-- alignment would need to bet against a richer alternative (the
@@ -30,7 +38,7 @@
-- @alpha = 1e-3@ against a stream of paired observations where @a@
-- runs systematically higher than @b@:
--
--- >>> let cfg = config 0.0 1.0 1.0e-3 Newton
+-- >>> let Right cfg = config 0.0 1.0 1.0e-3 Newton
-- >>> let ps = take 1000 (cycle [(1, 0), (1, 0), (0, 0), (1, 1)])
-- >>> decide cfg (foldl' (update cfg) (initial cfg) ps)
-- Reject
@@ -40,6 +48,7 @@ module Numeric.Eproc.Paired (
Config
, State
, Verdict(..)
+ , ConfigError(..)
-- * Bettor strategies
, Bettor(..)
@@ -58,7 +67,7 @@ module Numeric.Eproc.Paired (
) where
import qualified Numeric.Eproc.Bounded as Bounded
-import Numeric.Eproc.Common (Bettor(..), Verdict(..))
+import Numeric.Eproc.Common (Bettor(..), Verdict(..), ConfigError(..))
-- types ----------------------------------------------------------------------
@@ -80,16 +89,20 @@ newtype State = State Bounded.State
-- on the differences, which lie in @[lo - hi, hi - lo]@ with null
-- mean @0@.
--
--- >>> let cfg = config 0.0 1.0 1.0e-3 Newton
+-- Returns 'Left' with a 'ConfigError' on inputs that would leave
+-- the mathematical regime: @lo >= hi@ or @alpha@ outside
+-- @(0, 1)@.
+--
+-- >>> let Right cfg = config 0.0 1.0 1.0e-3 Newton
config
:: Double -- ^ sample lower bound @lo@
-> Double -- ^ sample upper bound @hi@
-> Double -- ^ significance level @alpha@
-> Bettor -- ^ bettor strategy
- -> Config
+ -> Either ConfigError Config
config !lo !hi !alpha b =
let !d = hi - lo
- in Config (Bounded.config 0 (negate d) d alpha b)
+ in fmap Config (Bounded.config 0 (negate d) d alpha b)
{-# INLINE config #-}
-- | The initial 'State' for a fresh streaming test.
@@ -106,6 +119,10 @@ initial (Config c) = State (Bounded.initial c)
-- Equivalent to feeding the difference @a - b@ into the underlying
-- bounded-mean test.
--
+-- /Precondition/: both @a@ and @b@ must lie in the @[lo, hi]@
+-- interval given to 'config'. The type-I error guarantee of the
+-- test depends on this; the function does not check.
+--
-- >>> let s1 = update cfg s0 (0.3, 0.7)
update :: Config -> State -> (Double, Double) -> State
update (Config c) (State s) (!a, !b) =
@@ -115,7 +132,7 @@ update (Config c) (State s) (!a, !b) =
-- | Compute the current 'Verdict' from the running 'State'.
--
-- 'Reject' iff either directional log-wealth of the underlying
--- bounded-mean test on the differences has crossed
+-- bounded-mean test on the differences has /ever/ crossed
-- @log(2 \/ alpha)@.
--
-- >>> decide cfg s0
@@ -126,8 +143,8 @@ decide (Config c) (State s) = Bounded.decide c s
-- inspection -----------------------------------------------------------------
--- | The current log-wealth of the underlying bounded-mean test on
--- the differences.
+-- | The supremum-so-far log-wealth of the underlying bounded-mean
+-- test on the differences.
--
-- >>> log_wealth s0
-- 0.0
diff --git a/test/Main.hs b/test/Main.hs
@@ -6,9 +6,11 @@ import Data.Bits
import Data.Word
import qualified Numeric.Eproc.Bernoulli as Bern
import qualified Numeric.Eproc.Bounded as Bounded
+import qualified Numeric.Eproc.Common as C
import qualified Numeric.Eproc.Paired as P
import Test.Tasty
import Test.Tasty.HUnit
+import qualified Test.Tasty.QuickCheck as QC
main :: IO ()
main = defaultMain $ testGroup "ppad-eproc" [
@@ -18,8 +20,16 @@ main = defaultMain $ testGroup "ppad-eproc" [
, two_sample_tests
, bernoulli_tests
, bettor_smoke_tests
+ , latched_rejection_tests
+ , config_validation_tests
+ , safety_property_tests
]
+-- partial helper: tests below hardcode valid configs.
+ok :: Either e a -> a
+ok (Right x) = x
+ok (Left _) = error "test: invalid config"
+
-- prng -----------------------------------------------------------------------
-- inline PCG-style PRNG, no external deps.
@@ -110,8 +120,8 @@ run_paired cfg pa pb budget g0 = go 0 g0 (P.initial cfg)
go !n !g !st
| n >= budget = (P.decide cfg st, n)
| otherwise = case P.decide cfg st of
- Bounded.Reject -> (Bounded.Reject, n)
- Bounded.Continue ->
+ P.Reject -> (P.Reject, n)
+ P.Continue ->
let (a, g1) = bernoulli pa g
(b, g2) = bernoulli pb g1
st' = P.update cfg st (a, b)
@@ -130,7 +140,7 @@ paired_avg_rate cfg pa pb budget trials seed =
rejects = length
[ () | g <- gens
, let (v, _) = run_paired cfg pa pb budget g
- , v == Bounded.Reject ]
+ , v == P.Reject ]
in fromIntegral rejects / fromIntegral trials
-- sanity ---------------------------------------------------------------------
@@ -139,12 +149,12 @@ paired_avg_rate cfg pa pb budget trials seed =
sanity_tests :: TestTree
sanity_tests = testGroup "sanity" [
testCase "degenerate input never rejects" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-6 Bounded.Newton
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-6 Bounded.Newton)
xs = replicate 5000 0.5
st = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs
Bounded.decide cfg st @?= Bounded.Continue
, testCase "two-sided thresholds applied symmetrically" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-6 Bounded.Newton
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-6 Bounded.Newton)
Bounded.decide cfg (Bounded.initial cfg) @?= Bounded.Continue
]
@@ -156,17 +166,17 @@ sanity_tests = testGroup "sanity" [
calibration_tests :: TestTree
calibration_tests = testGroup "null calibration" [
testCase "Newton, Bernoulli(0.5), m=0.5, alpha=0.05" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 0.05 Bounded.Newton
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 0.05 Bounded.Newton)
rate = rejection_rate cfg 0.5 2000 200 12345
- -- expected rate <= 0.05; allow up to 0.10 slack for sampling
- -- variability over 200 trials.
+ -- expected rate <= 0.05; allow up to ~0.08 slack for sampling
+ -- variability over 200 trials (sigma ~ 0.015).
assertBool ("FPR " ++ show rate ++ " exceeded slack") $
- rate <= 0.10
+ rate <= 0.08
, testCase "Adaptive, Bernoulli(0.5), m=0.5, alpha=0.05" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 0.05 Bounded.Adaptive
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 0.05 Bounded.Adaptive)
rate = rejection_rate cfg 0.5 2000 200 67890
assertBool ("FPR " ++ show rate ++ " exceeded slack") $
- rate <= 0.10
+ rate <= 0.08
]
-- power ----------------------------------------------------------------------
@@ -175,12 +185,12 @@ calibration_tests = testGroup "null calibration" [
power_tests :: TestTree
power_tests = testGroup "power" [
testCase "Newton detects Bernoulli(0.7) vs m=0.5" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
rate = rejection_rate cfg 0.7 5000 100 11111
assertBool ("power " ++ show rate ++ " too low") $
rate >= 0.95
, testCase "Adaptive detects Bernoulli(0.7) vs m=0.5" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive)
rate = rejection_rate cfg 0.7 5000 100 22222
assertBool ("power " ++ show rate ++ " too low") $
rate >= 0.95
@@ -191,11 +201,11 @@ power_tests = testGroup "power" [
two_sample_tests :: TestTree
two_sample_tests = testGroup "two-sample" [
testCase "identical distributions don't reject" $ do
- let cfg = P.config 0.0 1.0 1.0e-3 Bounded.Newton
+ let cfg = ok (P.config 0.0 1.0 1.0e-3 Bounded.Newton)
rate = paired_avg_rate cfg 0.5 0.5 2000 100 33333
assertBool ("FPR " ++ show rate) $ rate <= 0.05
, testCase "different distributions reject" $ do
- let cfg = P.config 0.0 1.0 1.0e-3 Bounded.Newton
+ let cfg = ok (P.config 0.0 1.0 1.0e-3 Bounded.Newton)
rate = paired_avg_rate cfg 0.3 0.7 5000 100 44444
assertBool ("power " ++ show rate) $ rate >= 0.95
]
@@ -238,27 +248,27 @@ bernoulli_rate cfg p budget trials seed =
bernoulli_tests :: TestTree
bernoulli_tests = testGroup "bernoulli" [
testCase "all-zero stream never rejects" $ do
- let cfg = Bern.config 1.0e-6 0.05 Bern.Newton
+ let cfg = ok (Bern.config 0.05 1.0e-6 Bern.Newton)
xs = replicate 5000 False
st = foldl' (Bern.update cfg) (Bern.initial cfg) xs
Bern.decide cfg st @?= Bern.Continue
, testCase "Newton FPR under H_0 (p = p_0 = 0.05)" $ do
- let cfg = Bern.config 0.05 0.05 Bern.Newton
+ let cfg = ok (Bern.config 0.05 0.05 Bern.Newton)
rate = bernoulli_rate cfg 0.05 2000 200 55555
assertBool ("FPR " ++ show rate ++ " exceeded slack") $
- rate <= 0.10
+ rate <= 0.08
, testCase "Adaptive FPR under H_0 (p = p_0 = 0.05)" $ do
- let cfg = Bern.config 0.05 0.05 Bern.Adaptive
+ let cfg = ok (Bern.config 0.05 0.05 Bern.Adaptive)
rate = bernoulli_rate cfg 0.05 2000 200 66666
assertBool ("FPR " ++ show rate ++ " exceeded slack") $
- rate <= 0.10
+ rate <= 0.08
, testCase "Newton detects p = 0.3 vs p_0 = 0.05" $ do
- let cfg = Bern.config 1.0e-3 0.05 Bern.Newton
+ let cfg = ok (Bern.config 0.05 1.0e-3 Bern.Newton)
rate = bernoulli_rate cfg 0.3 5000 100 77777
assertBool ("power " ++ show rate ++ " too low") $
rate >= 0.95
, testCase "Adaptive detects p = 0.3 vs p_0 = 0.05" $ do
- let cfg = Bern.config 1.0e-3 0.05 Bern.Adaptive
+ let cfg = ok (Bern.config 0.05 1.0e-3 Bern.Adaptive)
rate = bernoulli_rate cfg 0.3 5000 100 88888
assertBool ("power " ++ show rate ++ " too low") $
rate >= 0.95
@@ -270,19 +280,185 @@ bernoulli_tests = testGroup "bernoulli" [
-- deterministic stream.
bettor_smoke_tests :: TestTree
bettor_smoke_tests = testGroup "bettor smoke" [
- testCase "fixed bettor runs without error" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5)
+ testCase "fixed bettor runs without error (bounded)" $ do
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 (Bounded.Fixed 0.5))
xs = take 100 (cycle [0.0, 1.0])
st = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs
assertBool "samples advanced" (Bounded.samples st == 100)
- , testCase "Newton bettor runs without error" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton
+ , testCase "Newton bettor runs without error (bounded)" $ do
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Newton)
xs = take 100 (cycle [0.0, 1.0])
st = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs
assertBool "samples advanced" (Bounded.samples st == 100)
- , testCase "Adaptive bettor runs without error" $ do
- let cfg = Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive
+ , testCase "Adaptive bettor runs without error (bounded)" $ do
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 Bounded.Adaptive)
xs = take 100 (cycle [0.0, 1.0])
st = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs
assertBool "samples advanced" (Bounded.samples st == 100)
+ , testCase "fixed bettor runs without error (bernoulli)" $ do
+ let cfg = ok (Bern.config 0.5 1.0e-3 (Bern.Fixed 0.5))
+ xs = take 100 (cycle [True, False])
+ st = foldl' (Bern.update cfg) (Bern.initial cfg) xs
+ assertBool "samples advanced" (Bern.samples st == 100)
+ , testCase "Newton bettor runs without error (bernoulli)" $ do
+ let cfg = ok (Bern.config 0.5 1.0e-3 Bern.Newton)
+ xs = take 100 (cycle [True, False])
+ st = foldl' (Bern.update cfg) (Bern.initial cfg) xs
+ assertBool "samples advanced" (Bern.samples st == 100)
+ , testCase "Adaptive bettor runs without error (bernoulli)" $ do
+ let cfg = ok (Bern.config 0.5 1.0e-3 Bern.Adaptive)
+ xs = take 100 (cycle [True, False])
+ st = foldl' (Bern.update cfg) (Bern.initial cfg) xs
+ assertBool "samples advanced" (Bern.samples st == 100)
+ ]
+
+-- latched rejection ----------------------------------------------------------
+
+-- once the wealth crosses threshold, subsequent observations driving the
+-- current wealth back below threshold must not unrejection the test.
+latched_rejection_tests :: TestTree
+latched_rejection_tests = testGroup "latched rejection" [
+ testCase "bounded: cross then drown stays rejected" $ do
+ -- alpha = 0.5 => threshold log(2/0.5) = log 4 ~ 1.386.
+ -- Fixed 1.0 with x=1 grows log_w_pos by log 1.5 ~ 0.405/step;
+ -- five 1s push it past threshold. Then forty 0s drop it well
+ -- below.
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 0.5 (Bounded.Fixed 1.0))
+ xs1 = replicate 5 1.0
+ xs2 = replicate 40 0.0
+ st1 = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs1
+ st2 = foldl' (Bounded.update cfg) st1 xs2
+ Bounded.decide cfg st1 @?= Bounded.Reject
+ Bounded.decide cfg st2 @?= Bounded.Reject
+ , testCase "bernoulli: cross then drown stays rejected" $ do
+ let cfg = ok (Bern.config 0.05 0.5 (Bern.Fixed 1.0))
+ xs1 = replicate 5 True
+ xs2 = replicate 200 False
+ st1 = foldl' (Bern.update cfg) (Bern.initial cfg) xs1
+ st2 = foldl' (Bern.update cfg) st1 xs2
+ Bern.decide cfg st1 @?= Bern.Reject
+ Bern.decide cfg st2 @?= Bern.Reject
+ ]
+
+-- config validation ----------------------------------------------------------
+
+config_validation_tests :: TestTree
+config_validation_tests = testGroup "config validation" [
+ testCase "Bounded: alpha <= 0 rejected" $
+ assertLeft (Bounded.config 0.5 0.0 1.0 0.0 Bounded.Newton)
+ , testCase "Bounded: alpha >= 1 rejected" $
+ assertLeft (Bounded.config 0.5 0.0 1.0 1.5 Bounded.Newton)
+ , testCase "Bounded: lo >= hi rejected" $
+ assertLeft (Bounded.config 0.5 1.0 0.0 0.01 Bounded.Newton)
+ , testCase "Bounded: m == lo rejected" $
+ assertLeft (Bounded.config 0.0 0.0 1.0 0.01 Bounded.Newton)
+ , testCase "Bounded: m == hi rejected" $
+ assertLeft (Bounded.config 1.0 0.0 1.0 0.01 Bounded.Newton)
+ , testCase "Bounded: m outside [lo, hi] rejected" $
+ assertLeft (Bounded.config 2.0 0.0 1.0 0.01 Bounded.Newton)
+ , testCase "Bernoulli: alpha <= 0 rejected" $
+ assertLeft (Bern.config 0.5 0.0 Bern.Newton)
+ , testCase "Bernoulli: alpha >= 1 rejected" $
+ assertLeft (Bern.config 0.5 1.0 Bern.Newton)
+ , testCase "Bernoulli: p0 == 0 rejected" $
+ assertLeft (Bern.config 0.0 0.05 Bern.Newton)
+ , testCase "Bernoulli: p0 == 1 rejected" $
+ assertLeft (Bern.config 1.0 0.05 Bern.Newton)
+ , testCase "Paired: alpha out of range rejected" $
+ assertLeft (P.config 0.0 1.0 0.0 Bounded.Newton)
+ , testCase "Paired: lo >= hi rejected" $
+ assertLeft (P.config 1.0 0.0 0.01 Bounded.Newton)
+ ]
+ where
+ assertLeft :: Either C.ConfigError a -> Assertion
+ assertLeft e = case e of
+ Left _ -> pure ()
+ Right _ -> assertFailure "expected Left"
+
+-- safety properties ----------------------------------------------------------
+
+unit_double :: QC.Gen Double
+unit_double = QC.choose (0, 1)
+
+arb_bettor :: QC.Gen C.Bettor
+arb_bettor = QC.oneof [
+ pure C.Adaptive
+ , pure C.Newton
+ , C.Fixed <$> QC.choose (-10, 10) -- intentionally include unsafe values
+ ]
+
+finite :: Double -> Bool
+finite x = not (isNaN x) && not (isInfinite x)
+
+monotone_reject_bounded :: [Bounded.Verdict] -> Bool
+monotone_reject_bounded [] = True
+monotone_reject_bounded (Bounded.Continue : rest) = monotone_reject_bounded rest
+monotone_reject_bounded (Bounded.Reject : rest) = all (== Bounded.Reject) rest
+
+monotone_reject_bern :: [Bern.Verdict] -> Bool
+monotone_reject_bern [] = True
+monotone_reject_bern (Bern.Continue : rest) = monotone_reject_bern rest
+monotone_reject_bern (Bern.Reject : rest) = all (== Bern.Reject) rest
+
+safety_property_tests :: TestTree
+safety_property_tests = testGroup "safety properties" [
+ QC.testProperty "Bounded: log_wealth finite after any admissible stream" $
+ QC.forAll arb_bettor $ \b ->
+ QC.forAll (QC.listOf unit_double) $ \xs ->
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 b)
+ st = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs
+ in finite (Bounded.log_wealth st)
+
+ , QC.testProperty "Bernoulli: log_wealth finite after any admissible stream" $
+ QC.forAll arb_bettor $ \b ->
+ QC.forAll QC.arbitrary $ \xs ->
+ let cfg = ok (Bern.config 0.05 1.0e-3 b)
+ st = foldl' (Bern.update cfg) (Bern.initial cfg) (xs :: [Bool])
+ in finite (Bern.log_wealth st)
+
+ , QC.testProperty "Bounded: Fixed with arbitrary lambda is safe" $
+ QC.forAll (QC.choose (-1000, 1000)) $ \lam ->
+ QC.forAll (QC.listOf unit_double) $ \xs ->
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 (C.Fixed lam))
+ st = foldl' (Bounded.update cfg) (Bounded.initial cfg) xs
+ in finite (Bounded.log_wealth st)
+
+ , QC.testProperty "Bernoulli: Fixed with arbitrary lambda is safe" $
+ QC.forAll (QC.choose (-1000, 1000)) $ \lam ->
+ QC.forAll QC.arbitrary $ \xs ->
+ let cfg = ok (Bern.config 0.05 1.0e-3 (C.Fixed lam))
+ st = foldl' (Bern.update cfg) (Bern.initial cfg) (xs :: [Bool])
+ in finite (Bern.log_wealth st)
+
+ , QC.testProperty "Bounded: log_wealth is monotone nondecreasing" $
+ QC.forAll arb_bettor $ \b ->
+ QC.forAll (QC.listOf unit_double) $ \xs ->
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 1.0e-3 b)
+ sts = scanl (Bounded.update cfg) (Bounded.initial cfg) xs
+ lws = map Bounded.log_wealth sts
+ in and (zipWith (<=) lws (drop 1 lws))
+
+ , QC.testProperty "Bernoulli: log_wealth is monotone nondecreasing" $
+ QC.forAll arb_bettor $ \b ->
+ QC.forAll QC.arbitrary $ \xs ->
+ let cfg = ok (Bern.config 0.05 1.0e-3 b)
+ sts = scanl (Bern.update cfg) (Bern.initial cfg) (xs :: [Bool])
+ lws = map Bern.log_wealth sts
+ in and (zipWith (<=) lws (drop 1 lws))
+
+ , QC.testProperty "Bounded: rejection is latched" $
+ QC.forAll arb_bettor $ \b ->
+ QC.forAll (QC.listOf unit_double) $ \xs ->
+ let cfg = ok (Bounded.config 0.5 0.0 1.0 0.5 b)
+ sts = scanl (Bounded.update cfg) (Bounded.initial cfg) xs
+ vs = map (Bounded.decide cfg) sts
+ in monotone_reject_bounded vs
+
+ , QC.testProperty "Bernoulli: rejection is latched" $
+ QC.forAll arb_bettor $ \b ->
+ QC.forAll QC.arbitrary $ \xs ->
+ let cfg = ok (Bern.config 0.5 0.5 b)
+ sts = scanl (Bern.update cfg) (Bern.initial cfg) (xs :: [Bool])
+ vs = map (Bern.decide cfg) sts
+ in monotone_reject_bern vs
]