commit f60ba344661745c8829cf46c4539694086c473c2
parent a230e6558be4050da6ecb1a799d9993a9c957152
Author: Jared Tobin <jared@jtobin.io>
Date: Tue, 2 Jun 2026 21:39:33 -0230
style: align library with ppad conventions
Restyle the library to match the patterns used by ppad-fixed,
ppad-script, ppad-bip32, etc.
* snake_case for top-level functions ('log_wealth') and record fields
('cfg_bettor', 'st_n', 'cfg_lam_max_pos', 'st_log_w_pos', ...).
* Aligned record-field syntax matching bip32's HDKey.
* Section dividers ('-- types --------', '-- streaming --------', ...)
organising each module.
* '>>>' Haddock examples on every public binding; expanded haddocks
on the Bettor constructors describing each strategy and where
lambda_max comes from.
Also drop the Statistics.EProcess umbrella module: its renamed
camelCase exports ('meanConfig', 'updateMean', ...) existed only to
disambiguate overlapping names from Mean/TwoSample/Bettor and fought
the snake_case convention. Users now import the sub-modules directly,
mirroring ppad-fixed's layout:
import qualified Statistics.EProcess.Bettor as B
import qualified Statistics.EProcess.Mean as M
Tests and benches updated to the new imports and snake_case
conventions; README's GHCi example rewritten accordingly.
All 11 tests still pass; bench allocation and timing unchanged.
Diffstat:
10 files changed, 482 insertions(+), 445 deletions(-)
diff --git a/README.md b/README.md
@@ -1,40 +1,125 @@
-# ppad-eproc
+# eproc
-Anytime-valid sequential testing for Haskell, via e-processes and the
-betting framework.
+[](https://hackage.haskell.org/package/ppad-eproc)
+
+[](https://docs.ppad.tech/eproc)
-Implements the bounded-mean and paired two-sample tests of Waudby-Smith
-& Ramdas (2024) using predictable-plug-in bettors (fixed-lambda,
-aGRAPA, ONS). Tests are valid under optional stopping: reject as soon
-as the wealth process exceeds `1/alpha`, with type-I error controlled
-at `alpha` regardless of when you stop.
+Anytime-valid sequential hypothesis testing for bounded random
+variables, via the e-process / betting framework of
+[Waudby-Smith & Ramdas (2024)][wsr24]. Bounded-mean and paired
+two-sample tests are valid under optional stopping: reject as soon as
+the wealth process exceeds `1/alpha`, with type-I error controlled at
+`alpha` regardless of when the user stops streaming samples.
-## Use
+## Usage
-```haskell
-import qualified Statistics.EProcess as E
+A sample GHCi session:
--- Test H0: E[X] = 0.5 against H1: E[X] != 0.5,
--- samples bounded in [0, 1], alpha = 1e-6.
-let cfg = E.meanConfig 0.5 0.0 1.0 1.0e-6 E.ons
- s0 = E.initMeanState cfg
+```
+ > -- import qualified
+ > import qualified Statistics.EProcess.Bettor as B
+ > import qualified Statistics.EProcess.Mean as M
+ >
+ > -- test H_0: E[X] = 0.5 for samples in [0, 1] at alpha = 1e-3,
+ > -- with the ONS bettor
+ > let cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
+ >
+ > -- streaming interface: 'initial' then fold observations through 'update'
+ > let s0 = M.initial cfg
+ > let xs = [1, 1, 0, 1, 1, 0, 1, 1, 1, 1] -- mean 0.8, drifts from H_0
+ > let s10 = foldl (M.update cfg) s0 xs
+ >
+ > -- inspect wealth and verdict at any point
+ > M.samples s10
+ 10
+ > M.log_wealth s10
+ 0.7182493502552663
+ > M.decide cfg s10
+ Continue
+ >
+ > -- with enough evidence the test rejects
+ > let s300 = foldl (M.update cfg) s0 (concat (replicate 30 xs))
+ > M.log_wealth s300
+ 53.092214534054165
+ > M.decide cfg s300
+ Reject
+```
+
+For the paired two-sample mean-equality test, see
+`Statistics.EProcess.TwoSample`.
+
+## Documentation
+
+Haddocks (API documentation, etc.) are hosted at
+[docs.ppad.tech/eproc](https://docs.ppad.tech/eproc).
+
+## Performance
--- Stream samples through the test:
-let s1 = E.updateMean cfg s0 x1
- s2 = E.updateMean cfg s1 x2
- ...
+The aim is best-in-class performance for pure, highly-auditable Haskell
+code.
-case E.decideMean cfg sN of
- E.Reject -> ... -- H0 falsified
- E.Continue -> ... -- more data needed
+Current benchmark figures on an M4 Silicon MacBook Air look like (use
+`cabal bench` to run the benchmark suite):
+
+```
+ benchmarking Mean.update (one step)/ons
+ time 13.05 ns (12.95 ns .. 13.17 ns)
+ 1.000 R² (0.999 R² .. 1.000 R²)
+ mean 13.03 ns (12.95 ns .. 13.15 ns)
+ std dev 314.0 ps (248.3 ps .. 422.3 ps)
+
+ benchmarking Mean.update (1000-sample fold)/fixed
+ time 4.840 μs (4.819 μs .. 4.867 μs)
+ 1.000 R² (1.000 R² .. 1.000 R²)
+ mean 4.828 μs (4.817 μs .. 4.847 μs)
+ std dev 44.90 ns (30.94 ns .. 61.54 ns)
+
+ benchmarking Mean.update (1000-sample fold)/agrapa
+ time 15.67 μs (15.66 μs .. 15.69 μs)
+ 1.000 R² (1.000 R² .. 1.000 R²)
+ mean 15.67 μs (15.65 μs .. 15.69 μs)
+ std dev 63.74 ns (55.65 ns .. 75.07 ns)
+
+ benchmarking Mean.update (1000-sample fold)/ons
+ time 14.43 μs (14.42 μs .. 14.44 μs)
+ 1.000 R² (1.000 R² .. 1.000 R²)
+ mean 14.43 μs (14.42 μs .. 14.44 μs)
+ std dev 46.74 ns (34.00 ns .. 64.63 ns)
```
-For paired two-sample testing, see `Statistics.EProcess.TwoSample`.
+The inner update loop is fully fused: the `fixed` bettor allocates
+nothing per step, and the `agrapa` and `ons` bettors allocate a small
+constant per-step state record.
+
+You should compile with the `llvm` flag for maximum performance.
+
+## Development
+
+You'll require [Nix][nixos] with [flake][flake] support enabled. Enter a
+development shell with:
+
+```
+$ nix develop
+```
+
+Then do e.g.:
+
+```
+$ cabal repl ppad-eproc
+```
+
+to get a REPL for the main library.
+
+## References
-## Background
+- Waudby-Smith & Ramdas (2024), "[Estimating means of bounded random
+ variables by betting][wsr24]." JRSS-B.
+- Ramdas, Grunwald, Vovk, Shafer (2023), "[Game-theoretic statistics
+ and safe anytime-valid inference][rgvs23]." Statistical Science.
+- Shafer (2021), "[Testing by betting][shafer21]." JRSS-A.
-- Waudby-Smith & Ramdas (2024), "Estimating means of bounded random
- variables by betting." JRSS-B.
-- Ramdas, Grunwald, Vovk, Shafer (2023), "Game-theoretic statistics
- and safe anytime-valid inference." Statistical Science.
-- Shafer (2021), "Testing by betting." JRSS-A.
+[nixos]: https://nixos.org/
+[flake]: https://nixos.org/manual/nix/unstable/command-ref/new-cli/nix3-flake.html
+[wsr24]: https://arxiv.org/abs/2010.09686
+[rgvs23]: https://arxiv.org/abs/2210.01948
+[shafer21]: https://arxiv.org/abs/1909.03807
diff --git a/bench/Main.hs b/bench/Main.hs
@@ -4,7 +4,7 @@
module Main where
import Control.DeepSeq
-import qualified Statistics.EProcess as E
+import qualified Statistics.EProcess.Bettor as B
import qualified Statistics.EProcess.Mean as M
import qualified Statistics.EProcess.TwoSample as TS
import Criterion.Main
@@ -26,22 +26,22 @@ main = defaultMain [
update :: Benchmark
update =
- let !cfgF = M.config 0.5 0.0 1.0 1.0e-3 (E.Fixed 0.5)
- !cfgA = M.config 0.5 0.0 1.0 1.0e-3 E.Agrapa
- !cfgO = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
- !stF = M.initial cfgF
- !stA = M.initial cfgA
- !stO = M.initial cfgO
- !x = 0.7
+ let !cfg_f = M.config 0.5 0.0 1.0 1.0e-3 (B.Fixed 0.5)
+ !cfg_a = M.config 0.5 0.0 1.0 1.0e-3 B.Agrapa
+ !cfg_o = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
+ !st_f = M.initial cfg_f
+ !st_a = M.initial cfg_a
+ !st_o = M.initial cfg_o
+ !x = 0.7
in bgroup "Mean.update (one step)" [
- bench "fixed" $ nf (M.update cfgF stF) x
- , bench "agrapa" $ nf (M.update cfgA stA) x
- , bench "ons" $ nf (M.update cfgO stO) x
+ bench "fixed" $ nf (M.update cfg_f st_f) x
+ , bench "agrapa" $ nf (M.update cfg_a st_a) x
+ , bench "ons" $ nf (M.update cfg_o st_o) x
]
decide :: Benchmark
decide =
- let !cfg = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
+ let !cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
!st = M.initial cfg
in bgroup "Mean.decide" [
bench "initial state" $ nf (M.decide cfg) st
@@ -49,26 +49,26 @@ decide =
stream :: Benchmark
stream =
- let !xs = force (take 1000 (cycle [0.3, 0.7]))
- !cfgF = M.config 0.5 0.0 1.0 1.0e-3 (E.Fixed 0.5)
- !cfgA = M.config 0.5 0.0 1.0 1.0e-3 E.Agrapa
- !cfgO = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
- runM cfg = foldl' (M.update cfg) (M.initial cfg)
+ let !xs = force (take 1000 (cycle [0.3, 0.7]))
+ !cfg_f = M.config 0.5 0.0 1.0 1.0e-3 (B.Fixed 0.5)
+ !cfg_a = M.config 0.5 0.0 1.0 1.0e-3 B.Agrapa
+ !cfg_o = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
+ run_m cfg = foldl' (M.update cfg) (M.initial cfg)
in bgroup "Mean.update (1000-sample fold)" [
- bench "fixed" $ nf (runM cfgF) xs
- , bench "agrapa" $ nf (runM cfgA) xs
- , bench "ons" $ nf (runM cfgO) xs
+ bench "fixed" $ nf (run_m cfg_f) xs
+ , bench "agrapa" $ nf (run_m cfg_a) xs
+ , bench "ons" $ nf (run_m cfg_o) xs
]
twosample :: Benchmark
twosample =
- let !ps = force (take 1000 (cycle [(0.3, 0.7), (0.7, 0.3)]))
- !cfgF = TS.config 0.0 1.0 1.0e-3 (E.Fixed 0.5)
- !cfgA = TS.config 0.0 1.0 1.0e-3 E.Agrapa
- !cfgO = TS.config 0.0 1.0 1.0e-3 E.Ons
- runT cfg = foldl' (TS.update cfg) (TS.initial cfg)
+ let !ps = force (take 1000 (cycle [(0.3, 0.7), (0.7, 0.3)]))
+ !cfg_f = TS.config 0.0 1.0 1.0e-3 (B.Fixed 0.5)
+ !cfg_a = TS.config 0.0 1.0 1.0e-3 B.Agrapa
+ !cfg_o = TS.config 0.0 1.0 1.0e-3 B.Ons
+ run_t cfg = foldl' (TS.update cfg) (TS.initial cfg)
in bgroup "TwoSample.update (1000-sample fold)" [
- bench "fixed" $ nf (runT cfgF) ps
- , bench "agrapa" $ nf (runT cfgA) ps
- , bench "ons" $ nf (runT cfgO) ps
+ bench "fixed" $ nf (run_t cfg_f) ps
+ , bench "agrapa" $ nf (run_t cfg_a) ps
+ , bench "ons" $ nf (run_t cfg_o) ps
]
diff --git a/bench/Weight.hs b/bench/Weight.hs
@@ -4,7 +4,7 @@
module Main where
import Control.DeepSeq
-import qualified Statistics.EProcess as E
+import qualified Statistics.EProcess.Bettor as B
import qualified Statistics.EProcess.Mean as M
import qualified Statistics.EProcess.TwoSample as TS
import Weigh
@@ -23,44 +23,44 @@ main = mainWith $ do
update :: Weigh ()
update =
- let !cfgF = M.config 0.5 0.0 1.0 1.0e-3 (E.Fixed 0.5)
- !cfgA = M.config 0.5 0.0 1.0 1.0e-3 E.Agrapa
- !cfgO = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
- !stF = M.initial cfgF
- !stA = M.initial cfgA
- !stO = M.initial cfgO
+ let !cfg_f = M.config 0.5 0.0 1.0 1.0e-3 (B.Fixed 0.5)
+ !cfg_a = M.config 0.5 0.0 1.0 1.0e-3 B.Agrapa
+ !cfg_o = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
+ !st_f = M.initial cfg_f
+ !st_a = M.initial cfg_a
+ !st_o = M.initial cfg_o
in wgroup "Mean.update (one step)" $ do
- func "fixed" (M.update cfgF stF) 0.7
- func "agrapa" (M.update cfgA stA) 0.7
- func "ons" (M.update cfgO stO) 0.7
+ func "fixed" (M.update cfg_f st_f) 0.7
+ func "agrapa" (M.update cfg_a st_a) 0.7
+ func "ons" (M.update cfg_o st_o) 0.7
decide :: Weigh ()
decide =
- let !cfg = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
+ let !cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
!st = M.initial cfg
in wgroup "Mean.decide" $ do
func "initial state" (M.decide cfg) st
stream :: Weigh ()
stream =
- let !xs = force (take 1000 (cycle [0.3, 0.7]))
- !cfgF = M.config 0.5 0.0 1.0 1.0e-3 (E.Fixed 0.5)
- !cfgA = M.config 0.5 0.0 1.0 1.0e-3 E.Agrapa
- !cfgO = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
- runM cfg = foldl' (M.update cfg) (M.initial cfg)
+ let !xs = force (take 1000 (cycle [0.3, 0.7]))
+ !cfg_f = M.config 0.5 0.0 1.0 1.0e-3 (B.Fixed 0.5)
+ !cfg_a = M.config 0.5 0.0 1.0 1.0e-3 B.Agrapa
+ !cfg_o = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
+ run_m cfg = foldl' (M.update cfg) (M.initial cfg)
in wgroup "Mean.update (1000-sample fold)" $ do
- func "fixed" (runM cfgF) xs
- func "agrapa" (runM cfgA) xs
- func "ons" (runM cfgO) xs
+ func "fixed" (run_m cfg_f) xs
+ func "agrapa" (run_m cfg_a) xs
+ func "ons" (run_m cfg_o) xs
twosample :: Weigh ()
twosample =
- let !ps = force (take 1000 (cycle [(0.3, 0.7), (0.7, 0.3)]))
- !cfgF = TS.config 0.0 1.0 1.0e-3 (E.Fixed 0.5)
- !cfgA = TS.config 0.0 1.0 1.0e-3 E.Agrapa
- !cfgO = TS.config 0.0 1.0 1.0e-3 E.Ons
- runT cfg = foldl' (TS.update cfg) (TS.initial cfg)
+ let !ps = force (take 1000 (cycle [(0.3, 0.7), (0.7, 0.3)]))
+ !cfg_f = TS.config 0.0 1.0 1.0e-3 (B.Fixed 0.5)
+ !cfg_a = TS.config 0.0 1.0 1.0e-3 B.Agrapa
+ !cfg_o = TS.config 0.0 1.0 1.0e-3 B.Ons
+ run_t cfg = foldl' (TS.update cfg) (TS.initial cfg)
in wgroup "TwoSample.update (1000-sample fold)" $ do
- func "fixed" (runT cfgF) ps
- func "agrapa" (runT cfgA) ps
- func "ons" (runT cfgO) ps
+ func "fixed" (run_t cfg_f) ps
+ func "agrapa" (run_t cfg_a) ps
+ func "ons" (run_t cfg_o) ps
diff --git a/flake.nix b/flake.nix
@@ -16,27 +16,24 @@
let
lib = "ppad-eproc";
- pkgs = import nixpkgs { inherit system; };
- hlib = pkgs.haskell.lib;
- llvm = pkgs.llvmPackages_19.llvm;
- clang = pkgs.llvmPackages_19.clang;
+ pkgs = import nixpkgs { inherit system; };
+ hlib = pkgs.haskell.lib;
hpkgs = pkgs.haskell.packages.ghc910.extend (new: old: {
- ${lib} = new.callCabal2nix lib ./. {};
+ ${lib} = old.callCabal2nixWithOptions lib ./. "--enable-profiling" {};
});
+ cabal = hpkgs.cabal-install;
cc = pkgs.stdenv.cc;
ghc = hpkgs.ghc;
- cabal = hpkgs.cabal-install;
+ llvm = pkgs.llvmPackages_19.llvm;
in
{
packages.default = hpkgs.${lib};
- packages.haddock = hpkgs.${lib}.doc;
-
devShells.default = hpkgs.shellFor {
packages = p: [
- p.${lib}
+ (hlib.doBenchmark p.${lib})
];
buildInputs = [
@@ -45,12 +42,14 @@
llvm
];
+ doBenchmark = true;
+
shellHook = ''
PS1="[${lib}] \w$ "
echo "entering ${system} shell, using"
+ echo "cabal: $(${cabal}/bin/cabal --version)"
echo "cc: $(${cc}/bin/cc --version)"
echo "ghc: $(${ghc}/bin/ghc --version)"
- echo "cabal: $(${cabal}/bin/cabal --version)"
echo "llc: $(${llvm}/bin/llc --version | head -2 | tail -1)"
'';
};
diff --git a/lib/Statistics/EProcess.hs b/lib/Statistics/EProcess.hs
@@ -1,118 +0,0 @@
-{-# OPTIONS_HADDOCK prune #-}
-
--- |
--- Module: Statistics.EProcess
--- Copyright: (c) 2026 Jared Tobin
--- License: MIT
--- Maintainer: Jared Tobin <jared@ppad.tech>
---
--- Anytime-valid sequential hypothesis testing for bounded random
--- variables, via the e-process / betting framework of Waudby-Smith
--- and Ramdas (2024).
---
--- A bettor places predictable wagers against the null; the wealth
--- process is a nonnegative supermartingale under @H_0@, and Ville's
--- inequality gives type-I error control at @alpha@ for the stopping
--- rule \"reject the first time wealth exceeds @1\/alpha@\" —
--- regardless of when the user stops streaming samples.
---
--- This module re-exports the primary API. For finer control, see:
---
--- * "Statistics.EProcess.Bettor" for bettor strategies.
---
--- * "Statistics.EProcess.Mean" for the one-sample bounded-mean
--- test.
---
--- * "Statistics.EProcess.TwoSample" for the paired two-sample
--- mean-equality test.
-
-module Statistics.EProcess (
- -- * Bettors
- Bettor(..)
-
- -- * Bounded-mean test
- --
- -- $mean
- , Mean.Verdict(..)
- , meanConfig
- , initMeanState
- , updateMean
- , decideMean
-
- -- * Paired two-sample test
- --
- -- $twosample
- , twoSampleConfig
- , initTwoSampleState
- , updateTwoSample
- , decideTwoSample
-
- -- * Inspection
- , logWealth
- , samples
- ) where
-
-import Statistics.EProcess.Bettor
-import qualified Statistics.EProcess.Mean as Mean
-import qualified Statistics.EProcess.TwoSample as TS
-
--- $mean
---
--- For samples in @[lo, hi]@, test @H_0: E[x] = m@ two-sidedly.
-
--- | See 'Mean.config'.
-meanConfig
- :: Double -- ^ null mean @m@
- -> Double -- ^ sample lower bound
- -> Double -- ^ sample upper bound
- -> Double -- ^ significance level @alpha@
- -> Bettor
- -> Mean.Config
-meanConfig = Mean.config
-
--- | See 'Mean.initial'.
-initMeanState :: Mean.Config -> Mean.State
-initMeanState = Mean.initial
-
--- | See 'Mean.update'.
-updateMean :: Mean.Config -> Mean.State -> Double -> Mean.State
-updateMean = Mean.update
-
--- | See 'Mean.decide'.
-decideMean :: Mean.Config -> Mean.State -> Mean.Verdict
-decideMean = Mean.decide
-
--- $twosample
---
--- For paired observations @(a, b)@ both in @[lo, hi]@, test @H_0:
--- E[a] = E[b]@ two-sidedly.
-
--- | See 'TS.config'.
-twoSampleConfig
- :: Double
- -> Double
- -> Double
- -> Bettor
- -> TS.Config
-twoSampleConfig = TS.config
-
--- | See 'TS.initial'.
-initTwoSampleState :: TS.Config -> TS.State
-initTwoSampleState = TS.initial
-
--- | See 'TS.update'.
-updateTwoSample
- :: TS.Config -> TS.State -> (Double, Double) -> TS.State
-updateTwoSample = TS.update
-
--- | See 'TS.decide'.
-decideTwoSample :: TS.Config -> TS.State -> TS.Verdict
-decideTwoSample = TS.decide
-
--- | Current log-wealth of a 'Mean.State'.
-logWealth :: Mean.State -> Double
-logWealth = Mean.logWealth
-
--- | Sample count consumed so far.
-samples :: Mean.State -> Int
-samples = Mean.samples
diff --git a/lib/Statistics/EProcess/Bettor.hs b/lib/Statistics/EProcess/Bettor.hs
@@ -7,41 +7,51 @@
-- License: MIT
-- Maintainer: Jared Tobin <jared@ppad.tech>
--
--- Bettor strategies for the e-process framework. A bettor describes
--- how, given the history of centred observations @z = x - m@, the
--- next predictable bet @lambda@ is chosen.
+-- Bettor strategies for the e-process framework.
--
--- The bet placed at step @t@ depends only on data observed through
--- step @t-1@; this predictability is what makes the resulting wealth
--- process a nonnegative supermartingale under the null hypothesis,
--- and hence anytime-valid via Ville's inequality.
+-- A bettor describes how, given the history of centred observations
+-- @z = x - m@, the next predictable bet @lambda@ is chosen. The bet
+-- placed at step @t@ depends only on data observed through step @t-1@;
+-- this predictability is what makes the resulting wealth process a
+-- nonnegative supermartingale under the null hypothesis, and hence
+-- anytime-valid via Ville's inequality.
module Statistics.EProcess.Bettor (
- -- * Bettor strategies
+ -- * Bettor strategies
Bettor(..)
) where
+-- bettor strategies ----------------------------------------------------------
+
-- | A predictable bettor.
--
--- * 'Fixed' always bets the supplied @lambda@; useful for
--- smoke-testing the framework and as a numerical baseline.
+-- For 'Agrapa' and 'Ons', the safe bet bound @lambda_max@ is derived
+-- from the sample bounds supplied to the surrounding test
+-- configuration (e.g. 'Statistics.EProcess.Mean.config').
+--
+-- * 'Fixed' always bets the supplied @lambda@; useful for smoke
+-- testing the framework and as a numerical baseline.
--
-- * 'Agrapa' is the aGRAPA (approximate growth-rate adaptive
--- predictable plug-in) bettor; tracks empirical mean and variance
+-- predictable plug-in) bettor. Tracks empirical mean and variance
-- of centred observations and bets the Kelly-optimal value given
-- the current point estimate, clipped to @[0, lambda_max]@.
--
--- * 'Ons' is the online Newton step bettor; maintains a running
+-- * 'Ons' is the online Newton step bettor. Maintains a running
-- sum of squared gradients of the per-step log-wealth loss and
--- updates @lambda@ by a Newton step at each observation. Achieves
+-- updates @lambda@ by a Newton step at each observation; achieves
-- logarithmic regret against the best constant bet in hindsight,
--- and is in practice the strongest of the three under most signal
--- regimes.
+-- and is in practice the strongest of the three bettors under most
+-- signal regimes.
--
--- For 'Agrapa' and 'Ons', @lambda_max@ is derived from the sample
--- bounds supplied to the surrounding test 'Statistics.EProcess.Mean.config'.
-data Bettor
- = Fixed {-# UNPACK #-} !Double
+-- >>> Fixed 0.5
+-- Fixed 0.5
+-- >>> Agrapa
+-- Agrapa
+-- >>> Ons
+-- Ons
+data Bettor =
+ Fixed {-# UNPACK #-} !Double
| Agrapa
| Ons
deriving (Eq, Show)
diff --git a/lib/Statistics/EProcess/Mean.hs b/lib/Statistics/EProcess/Mean.hs
@@ -13,43 +13,47 @@
--
-- For samples @x_t@ in @[lo, hi]@, tests @H_0: E[x] = m@ against
-- @H_1: E[x] /= m@. Runs two e-processes simultaneously (one per
--- direction) and combines them by Bonferroni: reject if either
--- side's wealth crosses @2 \/ alpha@.
+-- direction) and combines them by Bonferroni: reject if either side's
+-- wealth crosses @2 \/ alpha@.
--
-- The test is anytime-valid: type-I error is controlled at @alpha@
-- regardless of when the user stops streaming samples.
module Statistics.EProcess.Mean (
- -- * Types
+ -- * Test configuration and state
Config
, State
, Verdict(..)
- -- * Construction
+ -- * Construction
, config
-
- -- * Streaming interface
, initial
+
+ -- * Streaming
, update
, decide
- -- * Inspection
- , logWealth
+ -- * Inspection
+ , log_wealth
, samples
) where
import GHC.Exts (Double(D#))
import Statistics.EProcess.Bettor
+-- types ----------------------------------------------------------------------
+
-- | Test outcome at the current sample count.
-data Verdict = Reject | Continue
+data Verdict =
+ Reject
+ | Continue
deriving (Eq, Show)
--- per-direction bettor state. one constructor per 'Bettor'
--- alternative; the constructor used in a given 'State' matches the
--- 'Bettor' chosen in the surrounding 'Config'.
-data BetState
- = SFixed
+-- 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
| SAgrapa
{-# UNPACK #-} !Double -- sum of z
{-# UNPACK #-} !Double -- sum of z^2
@@ -58,26 +62,28 @@ data BetState
{-# UNPACK #-} !Double -- lambda
{-# UNPACK #-} !Double -- acc (sum of squared gradients)
--- | Test configuration. Constructed by 'config'.
-data Config = Config
- { cfgBettor :: !Bettor
- , cfgLamMaxPos :: {-# UNPACK #-} !Double
- , cfgLamMaxNeg :: {-# UNPACK #-} !Double
- , cfgNullMean :: {-# UNPACK #-} !Double
- , cfgAlpha :: {-# UNPACK #-} !Double
- , cfgLogThresh :: {-# UNPACK #-} !Double
+-- | Bounded-mean test configuration. Build with 'config'.
+data Config = Config {
+ cfg_bettor :: !Bettor
+ , cfg_lam_max_pos :: {-# UNPACK #-} !Double
+ , cfg_lam_max_neg :: {-# UNPACK #-} !Double
+ , cfg_null_mean :: {-# UNPACK #-} !Double
+ , cfg_alpha :: {-# UNPACK #-} !Double
+ , cfg_log_thresh :: {-# UNPACK #-} !Double
}
--- | Test state. Two log-wealth processes (one per direction) and
--- per-direction bettor state.
-data State = State
- { stN :: {-# UNPACK #-} !Int
- , stLogWPos :: {-# UNPACK #-} !Double
- , stLogWNeg :: {-# UNPACK #-} !Double
- , stBetPos :: !BetState
- , stBetNeg :: !BetState
+-- | Streaming test state. Construct with 'initial' and fold
+-- observations through 'update'.
+data State = State {
+ st_n :: {-# UNPACK #-} !Int
+ , st_log_w_pos :: {-# UNPACK #-} !Double
+ , st_log_w_neg :: {-# UNPACK #-} !Double
+ , st_bet_pos :: !BetState
+ , st_bet_neg :: !BetState
}
+-- 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'
@@ -88,55 +94,18 @@ tiny :: Double
tiny = D# 1.0e-300##
{-# INLINE tiny #-}
--- | Build a test configuration.
---
--- >>> import qualified Statistics.EProcess.Bettor as B
--- >>> let cfg = config 0.5 0.0 1.0 1.0e-6 B.Ons
-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
- { cfgBettor = b
- , cfgLamMaxPos = 0.5 / (m - lo)
- , cfgLamMaxNeg = 0.5 / (hi - m)
- , cfgNullMean = m
- , cfgAlpha = alpha
- , cfgLogThresh = log (2 / alpha)
- }
--- NB. the lambda_max values are half the geometric ceiling; the 1/2
--- margin keeps the wealth factor bounded away from zero at the
--- boundary, which is the WSR safety recommendation.
-{-# INLINE config #-}
-
-- per-bettor initial state.
-initBet :: Bettor -> BetState
-initBet b = case b of
+init_bet :: Bettor -> BetState
+init_bet b = case b of
Fixed _ -> SFixed
Agrapa -> SAgrapa 0 0 0
Ons -> SOns 0 1.0e-6
-{-# INLINE initBet #-}
+{-# INLINE init_bet #-}
--- | Initial state for streaming.
-initial :: Config -> State
-initial Config{..} =
- let !s0 = initBet cfgBettor
- in State
- { stN = 0
- , stLogWPos = 0
- , stLogWNeg = 0
- , stBetPos = s0
- , stBetNeg = s0
- }
-{-# INLINE initial #-}
-
--- compute the next bet @lambda@ from the bettor and its current
--- state. @lamMax@ is the direction-specific safety bound.
-betLambda :: Bettor -> Double -> BetState -> Double
-betLambda b !lamMax !s = case b of
+-- compute the next bet 'lambda' from the bettor and its current
+-- state; 'lam_max' is the direction-specific safety bound.
+bet_lambda :: Bettor -> Double -> BetState -> Double
+bet_lambda b !lam_max !s = case b of
Fixed lam -> lam
Agrapa -> case s of
SAgrapa !sm !sm2 !n
@@ -148,16 +117,16 @@ betLambda b !lamMax !s = case b of
!var = max 0 (sm2 / nd - mu2)
!den = var + mu2
!raw = if den == 0 then 0 else mu / den
- in max 0 (min lamMax raw)
+ in max 0 (min lam_max raw)
_ -> 0
Ons -> case s of
SOns !lam _ -> lam
_ -> 0
-{-# INLINE betLambda #-}
+{-# INLINE bet_lambda #-}
--- update bettor state with newly observed centred value @z@.
-stepBet :: Bettor -> Double -> BetState -> Double -> BetState
-stepBet b !lamMax !s !z = case b of
+-- update bettor state with newly observed centred value 'z'.
+step_bet :: Bettor -> Double -> BetState -> Double -> BetState
+step_bet b !lam_max !s !z = case b of
Fixed _ -> SFixed
Agrapa -> case s of
SAgrapa !sm !sm2 !n -> SAgrapa (sm + z) (sm2 + z * z) (n + 1)
@@ -168,43 +137,100 @@ stepBet b !lamMax !s !z = case b of
!g = if denom == 0 then 0 else negate z / denom
!acc' = acc + g * g
!lam' = lam - g / acc'
- !clp = max 0 (min lamMax lam')
+ !clp = max 0 (min lam_max lam')
in SOns clp acc'
_ -> SOns 0 1.0e-6
-{-# INLINE stepBet #-}
+{-# INLINE step_bet #-}
+
+-- construction ---------------------------------------------------------------
+
+-- | Build a 'Config' for the bounded-mean test.
+--
+-- >>> import qualified Statistics.EProcess.Bettor as B
+-- >>> let cfg = config 0.5 0.0 1.0 1.0e-3 B.Ons
+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)
+ }
+-- NB. lambda_max values are half the geometric ceiling; the 1/2 margin
+-- keeps the wealth factor bounded away from zero at the boundary,
+-- which is the WSR safety recommendation.
+{-# INLINE config #-}
--- | Fold one observation into the state.
+-- | The initial 'State' for a fresh streaming test.
+--
+-- >>> let s0 = initial cfg
+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
+ }
+{-# INLINE initial #-}
+
+-- streaming ------------------------------------------------------------------
+
+-- | Fold one observation into the running 'State'.
+--
+-- >>> let s1 = update cfg s0 0.7
update :: Config -> State -> Double -> State
update Config{..} State{..} !x =
- let !z = x - cfgNullMean
- !lamP = betLambda cfgBettor cfgLamMaxPos stBetPos
- !lamN = betLambda cfgBettor cfgLamMaxNeg stBetNeg
- !facP = 1 + lamP * z
- !facN = 1 - lamN * z
- !logWP' = stLogWPos + log (max tiny facP)
- !logWN' = stLogWNeg + log (max tiny facN)
- !sP' = stepBet cfgBettor cfgLamMaxPos stBetPos z
- !sN' = stepBet cfgBettor cfgLamMaxNeg stBetNeg (negate z)
- in State (stN + 1) logWP' logWN' sP' sN'
+ let !z = x - cfg_null_mean
+ !lam_p = bet_lambda cfg_bettor cfg_lam_max_pos st_bet_pos
+ !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)
+ !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
{-# INLINE update #-}
--- | Decide based on current wealth.
+-- | Compute the current 'Verdict' from the running 'State'.
--
-- 'Reject' iff either directional log-wealth has crossed the
-- Bonferroni-adjusted threshold @log(2 \/ alpha)@.
+--
+-- >>> decide cfg s0
+-- Continue
decide :: Config -> State -> Verdict
decide Config{..} State{..}
- | stLogWPos >= cfgLogThresh = Reject
- | stLogWNeg >= cfgLogThresh = Reject
- | otherwise = Continue
+ | st_log_w_pos >= cfg_log_thresh = Reject
+ | st_log_w_neg >= cfg_log_thresh = Reject
+ | otherwise = Continue
{-# INLINE decide #-}
--- | Current log-wealth (the larger of the two directional processes).
-logWealth :: State -> Double
-logWealth State{..} = max stLogWPos stLogWNeg
-{-# INLINE logWealth #-}
+-- inspection -----------------------------------------------------------------
--- | Sample count consumed so far.
+-- | The current log-wealth, taken as the maximum of the two
+-- directional processes.
+--
+-- >>> log_wealth s0
+-- 0.0
+log_wealth :: State -> Double
+log_wealth State{..} = max st_log_w_pos st_log_w_neg
+{-# INLINE log_wealth #-}
+
+-- | The number of samples consumed so far.
+--
+-- >>> samples s0
+-- 0
samples :: State -> Int
-samples = stN
+samples = st_n
{-# INLINE samples #-}
diff --git a/lib/Statistics/EProcess/TwoSample.hs b/lib/Statistics/EProcess/TwoSample.hs
@@ -10,26 +10,26 @@
-- 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]@
--- by running the bounded-mean test on the differences @d_t = a_t -
--- b_t@ with null mean 0.
+-- @[lo, hi]@, tests @H_0: E[a] = E[b]@ against @H_1: E[a] /= E[b]@ by
+-- running the bounded-mean test on the differences @d_t = a_t - b_t@
+-- with null mean 0.
module Statistics.EProcess.TwoSample (
- -- * Types
+ -- * Test configuration and state
Config
, State
, Verdict(..)
- -- * Construction
+ -- * Construction
, config
-
- -- * Streaming interface
, initial
+
+ -- * Streaming
, update
, decide
- -- * Inspection
- , logWealth
+ -- * Inspection
+ , log_wealth
, samples
) where
@@ -37,19 +37,24 @@ import qualified Statistics.EProcess.Mean as M
import Statistics.EProcess.Mean (Verdict(..))
import Statistics.EProcess.Bettor (Bettor)
--- | Test configuration.
+-- types ----------------------------------------------------------------------
+
+-- | Paired two-sample test configuration. Build with 'config'.
newtype Config = Config M.Config
--- | Test state.
+-- | Streaming paired two-sample test state. Construct with 'initial'
+-- and fold observations through 'update'.
newtype State = State M.State
--- | Build a paired two-sample test configuration.
+-- construction ---------------------------------------------------------------
+
+-- | Build a 'Config' for the paired two-sample test.
--
-- Bounds @lo@ and @hi@ are the (shared) bounds on the individual
-- samples; differences then lie in @[lo - hi, hi - lo]@.
--
-- >>> import qualified Statistics.EProcess.Bettor as B
--- >>> let cfg = config 0.0 1.0 1.0e-6 B.Ons
+-- >>> let cfg = config 0.0 1.0 1.0e-3 B.Ons
config
:: Double -- ^ sample lower bound @lo@
-> Double -- ^ sample upper bound @hi@
@@ -61,28 +66,45 @@ config !lo !hi !alpha b =
in Config (M.config 0 (negate d) d alpha b)
{-# INLINE config #-}
--- | Initial state for streaming.
+-- | The initial 'State' for a fresh streaming test.
+--
+-- >>> let s0 = initial cfg
initial :: Config -> State
initial (Config c) = State (M.initial c)
{-# INLINE initial #-}
--- | Fold one paired observation @(a, b)@ into the state.
+-- streaming ------------------------------------------------------------------
+
+-- | Fold one paired observation @(a, b)@ into the running 'State'.
+--
+-- >>> let s1 = update cfg s0 (0.3, 0.7)
update :: Config -> State -> (Double, Double) -> State
update (Config c) (State s) (!a, !b) =
State (M.update c s (a - b))
{-# INLINE update #-}
--- | Decide based on current wealth.
+-- | Compute the current 'Verdict' from the running 'State'.
+--
+-- >>> decide cfg s0
+-- Continue
decide :: Config -> State -> Verdict
decide (Config c) (State s) = M.decide c s
{-# INLINE decide #-}
--- | Current log-wealth.
-logWealth :: State -> Double
-logWealth (State s) = M.logWealth s
-{-# INLINE logWealth #-}
+-- inspection -----------------------------------------------------------------
--- | Sample count consumed so far.
+-- | The current log-wealth.
+--
+-- >>> log_wealth s0
+-- 0.0
+log_wealth :: State -> Double
+log_wealth (State s) = M.log_wealth s
+{-# INLINE log_wealth #-}
+
+-- | The number of samples consumed so far.
+--
+-- >>> samples s0
+-- 0
samples :: State -> Int
samples (State s) = M.samples s
{-# INLINE samples #-}
diff --git a/ppad-eproc.cabal b/ppad-eproc.cabal
@@ -1,7 +1,7 @@
cabal-version: 3.0
name: ppad-eproc
version: 0.1.0
-synopsis: Anytime-valid sequential testing via e-processes
+synopsis: Anytime-valid sequential testing via e-processes.
license: MIT
license-file: LICENSE
author: Jared Tobin
@@ -34,7 +34,6 @@ library
if flag(llvm)
ghc-options: -fllvm -O2
exposed-modules:
- Statistics.EProcess
Statistics.EProcess.Bettor
Statistics.EProcess.Mean
Statistics.EProcess.TwoSample
diff --git a/test/Main.hs b/test/Main.hs
@@ -3,9 +3,8 @@
module Main where
import Data.Bits
-import Data.List (foldl')
import Data.Word
-import qualified Statistics.EProcess as E
+import qualified Statistics.EProcess.Bettor as B
import qualified Statistics.EProcess.Mean as M
import qualified Statistics.EProcess.TwoSample as TS
import Test.Tasty
@@ -13,47 +12,54 @@ import Test.Tasty.HUnit
main :: IO ()
main = defaultMain $ testGroup "ppad-eproc" [
- sanityTests
- , calibrationTests
- , powerTests
- , twoSampleTests
- , bettorSmokeTests
+ sanity_tests
+ , calibration_tests
+ , power_tests
+ , two_sample_tests
+ , bettor_smoke_tests
]
+-- prng -----------------------------------------------------------------------
+
-- inline PCG-style PRNG, no external deps.
newtype Gen = Gen Word64
-mkGen :: Word64 -> Gen
-mkGen = Gen
+mk_gen :: Word64 -> Gen
+mk_gen = Gen
-stepGen :: Gen -> (Word64, Gen)
-stepGen (Gen s) =
+step_gen :: Gen -> (Word64, Gen)
+step_gen (Gen s) =
let !s' = s * 6364136223846793005 + 1442695040888963407
in (s', Gen s')
-nextDouble :: Gen -> (Double, Gen)
-nextDouble g =
- let (w, g') = stepGen g
+next_double :: Gen -> (Double, Gen)
+next_double g =
+ let (w, g') = step_gen g
!x = fromIntegral (w `shiftR` 11 .&. 0x1FFFFFFFFFFFFF) /
9007199254740992
in (x, g')
bernoulli :: Double -> Gen -> (Double, Gen)
bernoulli !p g =
- let (u, g') = nextDouble g
+ let (u, g') = next_double g
in (if u < p then 1.0 else 0.0, g')
+gen_seq :: Gen -> [Gen]
+gen_seq g = let (_, g') = step_gen g in g : gen_seq g'
+
+-- harness --------------------------------------------------------------------
+
-- run a sequential mean test on a stream of n bernoulli(p) samples,
-- with the early-stopping rule built in. returns (verdict, samples
-- consumed).
-runMeanBernoulli
+run_mean_bernoulli
:: M.Config
-> Double -- ^ p
-> Int -- ^ budget
-> Gen
-> (M.Verdict, Int)
-runMeanBernoulli cfg p budget g0 = go 0 g0 (M.initial cfg)
+run_mean_bernoulli cfg p budget g0 = go 0 g0 (M.initial cfg)
where
go !n !g !st
| n >= budget = (M.decide cfg st, n)
@@ -65,133 +71,141 @@ runMeanBernoulli cfg p budget g0 = go 0 g0 (M.initial cfg)
in go (n + 1) g' st'
-- fraction of trials that rejected.
-rejectionRate
+rejection_rate
:: M.Config
-> Double -- ^ true bernoulli p
-> Int -- ^ budget per trial
-> Int -- ^ number of trials
-> Word64 -- ^ seed
-> Double
-rejectionRate cfg p budget trials seed =
- let gens = take trials (genSeq (mkGen seed))
+rejection_rate cfg p budget trials seed =
+ let gens = take trials (gen_seq (mk_gen seed))
rejects = length
[ () | g <- gens
- , let (v, _) = runMeanBernoulli cfg p budget g
+ , let (v, _) = run_mean_bernoulli cfg p budget g
, v == M.Reject ]
in fromIntegral rejects / fromIntegral trials
-genSeq :: Gen -> [Gen]
-genSeq g = let (_, g') = stepGen g in g : genSeq g'
+run_ts_paired
+ :: TS.Config
+ -> Double
+ -> Double -- ^ p for A and B
+ -> Int
+ -> Gen
+ -> (TS.Verdict, Int)
+run_ts_paired cfg pa pb budget g0 = go 0 g0 (TS.initial cfg)
+ where
+ go !n !g !st
+ | n >= budget = (TS.decide cfg st, n)
+ | otherwise = case TS.decide cfg st of
+ M.Reject -> (M.Reject, n)
+ M.Continue ->
+ let (a, g1) = bernoulli pa g
+ (b, g2) = bernoulli pb g1
+ st' = TS.update cfg st (a, b)
+ in go (n + 1) g2 st'
+
+ts_avg_rate
+ :: TS.Config
+ -> Double
+ -> Double
+ -> Int
+ -> Int
+ -> Word64
+ -> Double
+ts_avg_rate cfg pa pb budget trials seed =
+ let gens = take trials (gen_seq (mk_gen seed))
+ rejects = length
+ [ () | g <- gens
+ , let (v, _) = run_ts_paired cfg pa pb budget g
+ , v == M.Reject ]
+ in fromIntegral rejects / fromIntegral trials
--- sanity: with all-zero deviations from the null mean, no rejection.
+-- sanity ---------------------------------------------------------------------
-sanityTests :: TestTree
-sanityTests = testGroup "sanity" [
+-- with all-zero deviations from the null mean, no rejection.
+sanity_tests :: TestTree
+sanity_tests = testGroup "sanity" [
testCase "degenerate input never rejects" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-6 E.Ons
+ let cfg = M.config 0.5 0.0 1.0 1.0e-6 B.Ons
xs = replicate 5000 0.5
st = foldl' (M.update cfg) (M.initial cfg) xs
M.decide cfg st @?= M.Continue
, testCase "two-sided thresholds applied symmetrically" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-6 E.Ons
+ let cfg = M.config 0.5 0.0 1.0 1.0e-6 B.Ons
M.decide cfg (M.initial cfg) @?= M.Continue
]
--- null calibration: under H_0, with optional stopping, the empirical
--- rejection rate should be bounded by alpha. ville's inequality is
--- typically conservative on bernoulli, so the slack is small.
+-- null calibration -----------------------------------------------------------
-calibrationTests :: TestTree
-calibrationTests = testGroup "null calibration" [
+-- under H_0, with optional stopping, the empirical rejection rate should be
+-- bounded by alpha. ville's inequality is typically conservative on bernoulli,
+-- so the slack is small.
+calibration_tests :: TestTree
+calibration_tests = testGroup "null calibration" [
testCase "ONS, Bernoulli(0.5), m=0.5, alpha=0.05" $ do
- let cfg = M.config 0.5 0.0 1.0 0.05 E.Ons
- rate = rejectionRate cfg 0.5 2000 200 12345
- -- expected rate ≤ 0.05; allow up to 0.10 slack for sampling
+ let cfg = M.config 0.5 0.0 1.0 0.05 B.Ons
+ 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.
assertBool ("FPR " ++ show rate ++ " exceeded slack") $
rate <= 0.10
, testCase "aGRAPA, Bernoulli(0.5), m=0.5, alpha=0.05" $ do
- let cfg = M.config 0.5 0.0 1.0 0.05 E.Agrapa
- rate = rejectionRate cfg 0.5 2000 200 67890
+ let cfg = M.config 0.5 0.0 1.0 0.05 B.Agrapa
+ rate = rejection_rate cfg 0.5 2000 200 67890
assertBool ("FPR " ++ show rate ++ " exceeded slack") $
rate <= 0.10
]
--- power: under a clear shift, all (or nearly all) trials reject
--- within budget.
+-- power ----------------------------------------------------------------------
-powerTests :: TestTree
-powerTests = testGroup "power" [
+-- under a clear shift, all (or nearly all) trials reject within budget.
+power_tests :: TestTree
+power_tests = testGroup "power" [
testCase "ONS detects Bernoulli(0.7) vs m=0.5" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
- rate = rejectionRate cfg 0.7 5000 100 11111
+ let cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
+ rate = rejection_rate cfg 0.7 5000 100 11111
assertBool ("power " ++ show rate ++ " too low") $
rate >= 0.95
, testCase "aGRAPA detects Bernoulli(0.7) vs m=0.5" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-3 E.Agrapa
- rate = rejectionRate cfg 0.7 5000 100 22222
+ let cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Agrapa
+ rate = rejection_rate cfg 0.7 5000 100 22222
assertBool ("power " ++ show rate ++ " too low") $
rate >= 0.95
]
--- two-sample paired test.
+-- two-sample paired test -----------------------------------------------------
-runTSPaired
- :: TS.Config
- -> Double
- -> Double -- ^ p for A and B
- -> Int
- -> Gen
- -> (TS.Verdict, Int)
-runTSPaired cfg pA pB budget g0 = go 0 g0 (TS.initial cfg)
- where
- go !n !g !st
- | n >= budget = (TS.decide cfg st, n)
- | otherwise = case TS.decide cfg st of
- M.Reject -> (M.Reject, n)
- M.Continue ->
- let (a, g1) = bernoulli pA g
- (b, g2) = bernoulli pB g1
- st' = TS.update cfg st (a, b)
- in go (n + 1) g2 st'
-
-twoSampleTests :: TestTree
-twoSampleTests = testGroup "two-sample" [
+two_sample_tests :: TestTree
+two_sample_tests = testGroup "two-sample" [
testCase "identical distributions don't reject" $ do
- let cfg = TS.config 0.0 1.0 1.0e-3 E.Ons
- rate = avgRate cfg 0.5 0.5 2000 100 33333
+ let cfg = TS.config 0.0 1.0 1.0e-3 B.Ons
+ rate = ts_avg_rate cfg 0.5 0.5 2000 100 33333
assertBool ("FPR " ++ show rate) $ rate <= 0.05
, testCase "different distributions reject" $ do
- let cfg = TS.config 0.0 1.0 1.0e-3 E.Ons
- rate = avgRate cfg 0.3 0.7 5000 100 44444
+ let cfg = TS.config 0.0 1.0 1.0e-3 B.Ons
+ rate = ts_avg_rate cfg 0.3 0.7 5000 100 44444
assertBool ("power " ++ show rate) $ rate >= 0.95
]
- where
- avgRate cfg pA pB budget trials seed =
- let gens = take trials (genSeq (mkGen seed))
- rejects = length
- [ () | g <- gens
- , let (v, _) = runTSPaired cfg pA pB budget g
- , v == M.Reject ]
- in fromIntegral rejects / fromIntegral trials
-
--- bettor smoke tests: each bettor produces a well-defined state and
--- decision when run on a small deterministic stream.
-
-bettorSmokeTests :: TestTree
-bettorSmokeTests = testGroup "bettor smoke" [
+
+-- bettor smoke tests ---------------------------------------------------------
+
+-- each bettor produces a well-defined state and decision when run on a small
+-- deterministic stream.
+bettor_smoke_tests :: TestTree
+bettor_smoke_tests = testGroup "bettor smoke" [
testCase "fixed bettor runs without error" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-3 (E.Fixed 0.5)
+ let cfg = M.config 0.5 0.0 1.0 1.0e-3 (B.Fixed 0.5)
xs = take 100 (cycle [0.0, 1.0])
st = foldl' (M.update cfg) (M.initial cfg) xs
assertBool "samples advanced" (M.samples st == 100)
, testCase "ONS bettor runs without error" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-3 E.Ons
+ let cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Ons
xs = take 100 (cycle [0.0, 1.0])
st = foldl' (M.update cfg) (M.initial cfg) xs
assertBool "samples advanced" (M.samples st == 100)
, testCase "aGRAPA bettor runs without error" $ do
- let cfg = M.config 0.5 0.0 1.0 1.0e-3 E.Agrapa
+ let cfg = M.config 0.5 0.0 1.0 1.0e-3 B.Agrapa
xs = take 100 (cycle [0.0, 1.0])
st = foldl' (M.update cfg) (M.initial cfg) xs
assertBool "samples advanced" (M.samples st == 100)
