online 0.2.0 → 0.2.1.0
raw patch · 11 files changed
+702/−374 lines, 11 filesdep +doctestdep +formattingdep +onlinenew-component:exe:online-bench
Dependencies added: doctest, formatting, online, optparse-generic, perf, scientific, tasty, text
Files
- bench/bench.hs +213/−0
- online.cabal +80/−78
- readme.md +34/−4
- src/Online.hs +8/−6
- src/Online/Averages.hs +168/−0
- src/Online/Medians.hs +127/−0
- src/Online/Quantiles.hs +51/−31
- src/Online/Stats.hs +0/−140
- src/Online/StatsL1.hs +0/−112
- stack.yaml +2/−3
- test/test.hs +19/−0
+ bench/bench.hs view
@@ -0,0 +1,213 @@+{-# LANGUAGE NoImplicitPrelude #-}+{-# LANGUAGE OverloadedStrings #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE DeriveGeneric #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE FlexibleContexts #-}+{-# OPTIONS_GHC -Wall #-}+{-# OPTIONS_GHC -fno-warn-type-defaults #-}++import Online.Averages+import Online.Medians+import Options.Generic+import NumHask.Prelude hiding ((%))+import Perf+import Formatting+import qualified Data.Text as Text+import qualified Control.Foldl as L+import Data.Scientific++data Opts = Opts+ { runs :: Maybe Int -- <?> "number of runs"+ , sumTo :: Maybe Int -- <?> "sum to this number"+ } deriving (Generic, Show)++instance ParseRecord Opts++expt' :: Int -> Format r (Scientific -> r)+expt' x = scifmt Exponent (Just x)++code :: [Text] -> Text+code cs = "\n```\n" <> Text.intercalate "\n" cs <> "\n```\n"++sumInt :: [Int] -> Int -> Int+sumInt xs n = foldl' (+) zero (take n xs)++sumDouble :: [Double] -> Int -> Double+sumDouble xs n = foldl' (+) zero (take n xs)++sumPoly :: (Enum b, MultiplicativeUnital b, Additive b) => [b] -> Int -> b+sumPoly xs n = foldl' (+) zero (take n xs)++sumSum :: [Double] -> Int -> Double+sumSum xs n = L.fold L.sum (take n xs)++sumInt' :: Int -> Int+sumInt' x = foldl' (+) zero [zero .. x]++sumDouble' :: Double -> Double+sumDouble' x = foldl' (+) zero [one .. x]++sumPoly' :: (Enum b, MultiplicativeUnital b, Additive b) => b -> b+sumPoly' x = foldl' (+) zero [one .. x]++avTestMain :: [Double] -> Int -> Double+avTestMain xs n = L.fold Main.av (take n xs)++av :: (Field a) => L.Fold a a+av = L.Fold step begin extract+ where+ begin = (zero, zero)+ step (s, c) a = (s + a, c + one)+ extract (s, c) = s / c+{-# INLINABLE av #-}++maTest :: [Double] -> Int -> Double+maTest xs n = L.fold (ma 0.99) (take n xs)++formatRun :: [Cycle] -> Text -> Text+formatRun cs label =+ sformat+ ((right 24 ' ' %. stext) % stext %+ (left 7 ' ' %. expt' 3) % " cycles")+ label+ (Text.intercalate " " $ sformat (left 7 ' ' %. expt' 3) <$>+ (\x -> scientific (fromIntegral x) 0) <$> take 5 cs)+ (fromFloatDigits $ percentile 0.4 cs)++formatRunHeader :: Text+formatRunHeader =+ sformat+ ((right 24 ' ' %. stext) %+ (left 7 ' ' %. stext) %+ (left 8 ' ' %. stext) %+ (left 8 ' ' %. stext) %+ (left 8 ' ' %. stext) %+ (left 8 ' ' %. stext) %+ (left 8 ' ' %. stext))+ "run"+ "first"+ "2nd"+ "3rd"+ "4th"+ "5th"+ "40th %"++run :: Functor f => Text -> f ([Cycle], b) -> f Text+run label t = (`formatRun` label) . fst <$> t++tick_Test :: FilePath -> IO ()+tick_Test f = do+ onetick <- tick_+ ticks' <- replicateM 10 tick_+ manyticks <- replicateM 1000000 tick_+ let avticks = average manyticks+ let qticks = deciles 10 manyticks+ let tick999 = percentile 0.999 manyticks+ let tick99999 = percentile 0.99999 manyticks+ let tick99 = percentile 0.99 manyticks+ let tick40 = percentile 0.4 manyticks+ writeFile f $+ code+ [ "one tick_: " <> show onetick <> " cycles"+ , "next 10: " <> show ticks'+ , "average over 1m: " <> sformat (fixed 2) avticks <> " cycles"+ , "99.999% perc: " <> sformat commas (floor tick99999 :: Integer)+ , "99.9% perc: " <> sformat (fixed 2) tick999+ , "99th perc: " <> sformat (fixed 2) tick99+ , "40th perc: " <> sformat (fixed 2) tick40+ , "[min, 10th, 20th, .. 90th, max]:"+ , mconcat (sformat (" " % expt' 4) . fromFloatDigits <$> qticks)+ ]++fMono :: Int -> Int+fMono x = foldl' (+) 0 [1 .. x]++fLambda :: Int -> Int+fLambda = \x -> foldl' (+) 0 [1 .. x]++fPoly :: (Enum b, Num b, Additive b) => b -> b+fPoly x = foldl' (+) 0 [1 .. x]++tickTest :: FilePath -> Int -> IO ()+tickTest f a' = do+ (t, resultPrime) <- tick fMono a'+ print resultPrime+ (t2,_) <- tick fMono a'+ writeFile f $+ code+ [ "sum to " <> show a'+ , "first measure: " <> show t <> " cycles"+ , "second measure: " <> show t2 <> " cycles"+ ]+++ticksTest :: FilePath -> Int -> Int -> IO ()+ticksTest f a' n = do+ -- | various versions of tick+ rpure <- run "ticks" $ ticks n fMono a'+ rpurePoly <- run "ticks (poly)" $ ticks n fPoly a'+ rpureLambda <- run "ticks (lambda)" $ ticks n fLambda a'+ rio <- run "ticksIO" $ ticksIO n (pure $ fMono a')+ rioPoly <- run "ticksIO (poly)" $ ticksIO n (pure $ fPoly a')+ rioLambda <- run "ticksIO (lambda)" $ ticksIO n (pure $ fLambda a')++ writeFile f $+ code [ "sum to " <> show a' <> " n = " <> show n+ , formatRunHeader+ , rpure+ , rpureLambda+ , rpurePoly+ , rio+ , rioLambda+ , rioPoly+ ]++main :: IO ()+main = do+ o :: Opts <- getRecord "online performance benchmarking"+ let n = fromMaybe 100 (runs o)+ let a = fromMaybe 1000 (sumTo o)+ let !a' = fromIntegral a :: Double+ let !xs = [zero .. a]+ let !xs' = fromIntegral <$> xs :: [Double]+ _ <- warmup 100++ rSumInt' <- run "sumInt [0..]" $ ticks n sumInt' a+ rSumDouble' <- run "sumDouble [0..]" $ ticks n sumDouble' a'+ rSumPoly' <- run "sumPoly [0..]" $ ticks n sumPoly' a'+ rSumInt <- run "sum Int" $ ticks n (sumInt xs) a+ rSumDouble <- run "sum Double" $ ticks n (sumDouble xs') a+ rSumPoly <- run "sum Poly" $ ticks n (sumPoly xs') a+ rSumSum <- run "fold sum" $ ticks n (sumSum xs') a+ rAvTestMain <- run "fold av" $ ticks n (avTestMain xs') a+ rMaTest <- run "fold ma" $ ticks n (maTest xs') a+ rStdTest <- run "fold std" $ ticks n (\x -> L.fold (std 0.99) $ take x xs') a++ rMaL1Test <- run "fold maL1" $ ticks n (\x -> L.fold (maL1 0 0.01 0.99) $ take x xs') a+ rabsmaL1Test <- run "fold absmaL1" $ ticks n (\x -> L.fold (absmaL1 0 0.01 0.99) $ take x xs') a++ writeFile "other/perf.md" $+ code+ [ "sum to " <> sformat commas a+ , formatRunHeader+ , rSumInt'+ , rSumDouble'+ , rSumPoly'+ , rSumInt+ , rSumDouble+ , rSumPoly+ , rSumSum+ , rAvTestMain+ , rMaTest+ , rStdTest+ , rMaL1Test+ , rabsmaL1Test+ ]++ tick_Test "other/tick_.md"+ tickTest "other/tick.md" a+ ticksTest "other/ticks.md" a n+ pure ()
online.cabal view
@@ -1,85 +1,87 @@-name: online-version: 0.2.0-synopsis:- online statistics-description:- transformation of statistics to online algorithms-homepage:- https://github.com/tonyday567/online-license:- BSD3-license-file:- LICENSE-author:- Tony Day-maintainer:- tonyday567@gmail.com-copyright:- Tony Day-category:- statistics -build-type:- Simple-cabal-version:- >=1.10+-- This file has been generated from package.yaml by hpack version 0.20.0.+--+-- see: https://github.com/sol/hpack+--+-- hash: 160bd8de21b13879993819afe4d9b4640f9183a6e979143c22b2aaa8e4510e31++name: online+version: 0.2.1.0+synopsis: online statistics+description: transformation of statistics to online algorithms+category: statistics+homepage: https://github.com/tonyday567/online#readme+bug-reports: https://github.com/tonyday567/online/issues+author: Tony Day+maintainer: tonyday567@gmail.com+copyright: Tony Day+license: BSD3+license-file: LICENSE+build-type: Simple+cabal-version: >= 1.10+ extra-source-files:- readme.md- stack.yaml+ readme.md+ stack.yaml +source-repository head+ type: git+ location: https://github.com/tonyday567/online+ library- default-language:- Haskell2010- ghc-options:- -Wall hs-source-dirs:- src+ src+ default-extensions: NegativeLiterals NoImplicitPrelude OverloadedStrings UnicodeSyntax+ ghc-options: -Wall+ build-depends:+ base >=4.7 && <5+ , foldl+ , numhask+ , protolude+ , tdigest+ , vector+ , vector-algorithms exposed-modules:- Online,- Online.Stats,- Online.StatsL1,- Online.Quantiles+ Online+ Online.Averages+ Online.Medians+ Online.Quantiles+ other-modules:+ Paths_online+ default-language: Haskell2010++executable online-bench+ main-is: bench.hs+ hs-source-dirs:+ bench+ default-extensions: NegativeLiterals NoImplicitPrelude OverloadedStrings UnicodeSyntax BangPatterns BinaryLiterals DeriveFoldable DeriveFunctor DeriveGeneric DeriveTraversable DisambiguateRecordFields EmptyCase FlexibleContexts FlexibleInstances FunctionalDependencies GADTSyntax InstanceSigs KindSignatures LambdaCase MonadComprehensions MultiParamTypeClasses MultiWayIf ParallelListComp PartialTypeSignatures PatternSynonyms RankNTypes RecordWildCards RecursiveDo ScopedTypeVariables TupleSections TypeFamilies TypeOperators+ ghc-options: -O2 -funbox-strict-fields -fforce-recomp -threaded -rtsopts -with-rtsopts=-N build-depends:- base >= 4.7 && < 5,- protolude,- foldl,- vector,- tdigest,- numhask,- vector-algorithms- default-extensions:- NoImplicitPrelude,- UnicodeSyntax,- BangPatterns,- BinaryLiterals,- DeriveFoldable,- DeriveFunctor,- DeriveGeneric,- DeriveTraversable,- DisambiguateRecordFields,- EmptyCase,- FlexibleContexts,- FlexibleInstances,- FunctionalDependencies,- GADTSyntax,- InstanceSigs,- KindSignatures,- LambdaCase,- MonadComprehensions,- MultiParamTypeClasses,- MultiWayIf,- NegativeLiterals,- OverloadedStrings,- ParallelListComp,- PartialTypeSignatures,- PatternSynonyms,- RankNTypes,- RecordWildCards,- RecursiveDo,- ScopedTypeVariables,- TupleSections,- TypeFamilies,- TypeOperators+ base >=4.7 && <5+ , foldl+ , formatting+ , numhask+ , online+ , optparse-generic+ , perf+ , protolude+ , scientific+ , text+ other-modules:+ Paths_online+ default-language: Haskell2010 -source-repository head- type: git- location: https://github.com/tonyday567/online+test-suite test+ type: exitcode-stdio-1.0+ main-is: test.hs+ hs-source-dirs:+ test+ default-extensions: NegativeLiterals NoImplicitPrelude OverloadedStrings UnicodeSyntax+ build-depends:+ base >=4.7 && <5+ , doctest+ , online+ , protolude+ , tasty+ other-modules:+ Paths_online+ default-language: Haskell2010
readme.md view
@@ -1,12 +1,17 @@ [online](https://github.com/tonyday567/online)-===+============================================== -[](https://travis-ci.org/tonyday567/online) [](https://hackage.haskell.org/package/online) [](http://stackage.org/lts/package/online) [](http://stackage.org/nightly/package/online)+[](https://travis-ci.org/tonyday567/online)+[](https://hackage.haskell.org/package/online)+[](http://stackage.org/lts/package/online)+[](http://stackage.org/nightly/package/online) -online turns a statistic (in haskell this can usually be thought of as a fold of a foldable) into an online algorithm.+online turns a statistic (in haskell this can usually be thought of as a+fold of a foldable) into an online algorithm. motivation-===+========== Imagine a data stream, like an ordered indexed container or a time-series of measurements. An exponential moving average can be@@ -48,3 +53,28 @@ where ma r = online identity (.*r) sqma r = online (**2) (.*r)++[perf](https://hackage.haskell.org/package/perf)+================================================++1 cycle = 0.4 nanoseconds.++ sum to 1,000+ run first 2nd 3rd 4th 5th 40th %+ sumInt [0..] 8.69e3 1.59e3 1.54e3 1.60e3 1.64e3 1.57e3 cycles+ sumDouble [0..] 3.00e5 2.82e5 3.27e5 2.99e5 2.61e5 9.38e4 cycles+ sumPoly [0..] 9.34e4 1.05e5 1.09e5 1.18e5 1.64e5 7.81e4 cycles+ sum Int 1.67e4 1.17e4 1.16e4 1.17e4 1.16e4 1.17e4 cycles+ sum Double 2.63e4 1.18e4 1.17e4 1.16e4 1.16e4 1.16e4 cycles+ sum Poly 1.28e4 1.18e4 1.17e4 1.17e4 1.18e4 1.17e4 cycles+ fold sum 1.23e4 1.19e4 1.18e4 1.17e4 1.18e4 1.18e4 cycles+ fold av 1.23e4 1.18e4 1.18e4 1.18e4 1.18e4 1.18e4 cycles+ fold ma 1.24e4 1.19e4 1.20e4 1.19e4 1.19e4 1.19e4 cycles+ fold std 2.04e5 1.15e5 6.14e5 1.16e5 1.16e5 1.12e5 cycles+ fold maL1 9.11e4 8.32e4 1.27e5 3.34e5 1.04e5 1.02e5 cycles+ fold absmaL1 6.82e4 6.71e4 6.64e4 6.66e4 2.97e5 6.65e4 cycles++recipe+======++ stack build --test --exec "$(stack path --local-install-root)/bin/online-bench" --exec "$(stack path --local-bin)/pandoc -f markdown -i other/header.md other/readme_.md other/footer.md -t html -o index.html --filter pandoc-include --mathjax" --exec "$(stack path --local-bin)/pandoc -f markdown -i other/readme_.md -t markdown -o readme.md --filter pandoc-include --mathjax" --file-watch
src/Online.hs view
@@ -1,8 +1,10 @@ -- | online library-module Online (- module X-) where+module Online+ ( module Online.Averages+ , module Online.Medians+ , module Online.Quantiles+ ) where -import Online.Stats as X-import Online.StatsL1 as X-import Online.Quantiles as X+import Online.Quantiles+import Online.Averages+import Online.Medians
+ src/Online/Averages.hs view
@@ -0,0 +1,168 @@+{-# LANGUAGE NoImplicitPrelude #-}++-- | online statistics based on a moving average+module Online.Averages+ ( Averager+ , online+ -- * online statistics+ , ma+ , absma+ , sqma+ , std+ , cov+ , corr+ , corrGauss+ , beta+ , alpha+ , autocorr+ , mconst+ ) where++import qualified Control.Foldl as L+import Control.Foldl (Fold(..))+import NumHask.Prelude++-- | Most common statistics are averages.+newtype Averager a b = Averager+ { _averager :: (a, b)+ }++instance (Monoid a, Monoid b) => Monoid (Averager a b) where+ mempty = Averager (mempty, mempty)+ mappend (Averager (s, c)) (Averager (s', c')) =+ Averager (mappend s s', mappend c c')++-- | online takes a function and turns it into a `Control.Foldl.Fold` where the step is an incremental update of the (isomorphic) statistic.+online :: (Field b) => (a -> b) -> (b -> b) -> Fold a b+online f g = Fold step begin extract+ where+ begin = Averager (zero, zero)+ step (Averager (s, c)) a = Averager (g $ s + f a, g $ c + one)+ extract (Averager (s, c)) = s / c+{-# INLINABLE online #-}++-- $setup+--+-- >>> :set -XNoImplicitPrelude+-- >>> import NumHask.Prelude+-- >>> import qualified Control.Foldl as L+-- >>> let n = 100+-- >>> let r = 0.9++-- | moving average with a decay rate+-- +-- so 'ma 1' is the simple average (no decay in the statistic), and 'ma 0.00001' is the last value (insta-decay)+--+-- >>> L.fold (ma 1) [0..100]+-- 50.0+--+-- >>> L.fold (ma 1e-12) [0..100] ≈ 100+-- True+--+-- >>> L.fold (ma 0.9) [0..100]+-- 91.00241448887785+--+ma :: (Field a) => a -> Fold a a+ma r = online identity (* r)+{-# INLINABLE ma #-}++-- | absolute average+absma :: (Field a, Signed a) => a -> Fold a a+absma r = online abs (* r)+{-# INLINABLE absma #-}++-- | average square+sqma :: (Field a) => a -> Fold a a+sqma r = online (\x -> x * x) (* r)+{-# INLINABLE sqma #-}++-- | standard deviation+--+-- The formulae for standard deviation, expressed in online terminology, highlights how this statistic is composed of averages:+--+-- > (\s ss -> sqrt (ss - s ** (one+one))) <$> ma r <*> sqma r+--+-- The average deviation of the numbers 1..1000 is about 1 / sqrt 12 * 1000 (see <<https://en.wikipedia.org/wiki/Uniform_distribution_(continuous)#Standard_uniform wiki>>)+--+-- >>> L.fold (std 1) [0..1000]+-- 288.9636655359978+--+-- The average deviation with a decay of 0.99+--+-- >>> L.fold (std 0.99) [0..1000]+-- 99.28328803164005+std :: (ExpField a) => a -> Fold a a+std r = (\s ss -> sqrt (ss - s ** (one+one))) <$> ma r <*> sqma r+{-# INLINABLE std #-}++-- | the covariance of a tuple+-- given an underlying central tendency fold+cov :: (Field a) => Fold a a -> Fold (a, a) a+cov m =+ (\xy x' y' -> xy - x' * y') <$> L.premap (uncurry (*)) m <*> L.premap fst m <*>+ L.premap snd m+{-# INLINABLE cov #-}++-- | correlation of a tuple, specialised to Guassian+corrGauss :: (ExpField a) => a -> Fold (a, a) a+corrGauss r =+ (\cov' stdx stdy -> cov' / (stdx * stdy)) <$> cov (ma r) <*>+ L.premap fst (std r) <*>+ L.premap snd (std r)+{-# INLINABLE corrGauss #-}++-- | a generalised version of correlation of a tuple+corr :: (Field a) => Fold a a -> Fold a a -> Fold (a, a) a+corr central deviation =+ (\cov' stdx stdy -> cov' / (stdx * stdy)) <$> cov central <*>+ L.premap fst deviation <*>+ L.premap snd deviation+{-# INLINABLE corr #-}++-- | the beta in a simple linear regression of a tuple+-- given an underlying central tendency fold+beta :: (Field a) => Fold a a -> Fold (a, a) a+beta m =+ (\xy x' y' x2 -> (xy - x' * y') / (x2 - x' * x')) <$> L.premap (uncurry (*)) m <*>+ L.premap fst m <*>+ L.premap snd m <*>+ L.premap (\(x, _) -> x * x) m+{-# INLINABLE beta #-}++-- | the alpha of a tuple+alpha :: (Field a) => Fold a a -> Fold (a, a) a+alpha m = (\y b x -> y - b * x) <$> L.premap fst m <*> beta m <*> L.premap snd m+{-# INLINABLE alpha #-}++{-| autocorrelation is a slippery concept. This method starts with the concept that there is an underlying random error process (e), and autocorrelation is a process on top of that ie for a one-step correlation relationship.++value@t = e@t + k * e@t-1++where k is the autocorrelation.++There are thus two online rates needed: one for the average being considered to be the dependent variable, and one for the online of the correlation calculation between the most recent value and the moving average. For example,++> L.fold (autocorr zero one)++would estimate the one-step autocorrelation relationship of the previous value and the current value over the entire sample set.++-}+autocorr :: (BoundedField a) => Fold a a -> Fold (a, a) a -> Fold a a+autocorr central corrf =+ case central of+ (Fold mStep mBegin mDone) ->+ case corrf of+ (Fold dStep dBegin dDone) ->+ let begin = (mBegin, dBegin)+ step (mAcc, dAcc) a =+ ( mStep mAcc a+ , if isNaN (mDone mAcc)+ then dAcc+ else dStep dAcc (mDone mAcc, a))+ done = dDone . snd+ in Fold step begin done+{-# INLINABLE autocorr #-}++-- | a constant fold+mconst :: (Field a) => a -> L.Fold a a+mconst a = L.Fold (\() _ -> ()) () (const a)
+ src/Online/Medians.hs view
@@ -0,0 +1,127 @@+module Online.Medians+ ( -- * convert a statistic to an online median stat equivalent to L1+ Medianer(..)+ , onlineL1+ , onlineL1'++ -- * online statistics+ , maL1+ , absmaL1+ , covL1+ , corrL1+ , betaL1+ , alphaL1+ , autocorrL1+ ) where++import qualified Control.Foldl as L+import Control.Foldl (Fold(..))+import Protolude++-- | A rough Median.+-- The average absolute value of the stat is used to callibrate estimate drift towards the median+data Medianer a b = Medianer+ { medAbsSum :: a+ , medCount :: b+ , medianEst :: a+ }++-- | onlineL1' takes a function and turns it into a `Control.Foldl.Fold` where the step is an incremental update of an (isomorphic) median statistic.+onlineL1' ::+ (Ord b, Fractional b) => b -> b -> (a -> b) -> (b -> b) -> Fold a (b, b)+onlineL1' i d f g = Fold step begin extract+ where+ begin = Medianer 0 0 0+ step (Medianer s c m) a =+ Medianer+ (g $ s + abs (f a))+ (g $ c + 1)+ ((1 - d) * (m + s' * i * s / c') + d * f a)+ where+ c' =+ if c == 0+ then 1+ else c+ s'+ | f a > m = 1+ | f a < m = -1+ | otherwise = 0+ extract (Medianer s c m) = (s / c, m)+{-# INLINABLE onlineL1' #-}++-- | onlineL1 takes a function and turns it into a `Control.Foldl.Fold` where the step is an incremental update of an (isomorphic) median statistic.+onlineL1 :: (Ord b, Fractional b) => b -> b -> (a -> b) -> (b -> b) -> Fold a b+onlineL1 i d f g = snd <$> onlineL1' i d f g+{-# INLINABLE onlineL1 #-}++-- $setup+--+-- >>> :set -XNoImplicitPrelude+-- >>> import NumHask.Prelude+-- >>> import qualified Control.Foldl as L+-- >>> let n = 100+-- >>> let inc = 0.1+-- >>> let d = 0+-- >>> let r = 0.9++-- | moving median+-- >>> L.fold (maL1 inc d r) [1..n]+-- 93.92822312742108+--+maL1 :: (Ord a, Fractional a) => a -> a -> a -> Fold a a+maL1 i d r = onlineL1 i d identity (* r)+{-# INLINABLE maL1 #-}++-- | moving absolute deviation+absmaL1 :: (Ord a, Fractional a) => a -> a -> a -> Fold a a+absmaL1 i d r = fst <$> onlineL1' i d identity (* r)+{-# INLINABLE absmaL1 #-}++-- | covariance of a tuple+covL1 :: (Ord a, Fractional a) => a -> a -> a -> Fold (a, a) a+covL1 i d r =+ (\xy xbar ybar -> xy - xbar * ybar) <$> onlineL1 i d (uncurry (*)) (* r) <*>+ onlineL1 i d fst (* r) <*>+ onlineL1 i d snd (* r)+{-# INLINABLE covL1 #-}++-- | correlation of a tuple+corrL1 :: (Ord a, Floating a) => a -> a -> a -> Fold (a, a) a+corrL1 i d r =+ (\cov' stdx stdy -> cov' / (stdx * stdy)) <$> covL1 i d r <*>+ L.premap fst (absmaL1 i d r) <*>+ L.premap snd (absmaL1 i d r)+{-# INLINABLE corrL1 #-}++-- | the beta in a simple linear regression of a tuple+betaL1 :: (Ord a, Floating a) => a -> a -> a -> Fold (a, a) a+betaL1 i d r =+ (\xy x' y' x2 -> (xy - x' * y') / (x2 - x' * x')) <$>+ L.premap (uncurry (*)) (maL1 i d r) <*>+ L.premap fst (maL1 i d r) <*>+ L.premap snd (maL1 i d r) <*>+ L.premap (\(x, _) -> x * x) (maL1 i d r)+{-# INLINABLE betaL1 #-}++-- | the alpha in a simple linear regression of `snd` on `fst`+alphaL1 :: (Ord a, Floating a) => a -> a -> a -> Fold (a, a) a+alphaL1 i d r =+ (\y b x -> y - b * x) <$> L.premap fst (maL1 i d r) <*> betaL1 i d r <*>+ L.premap snd (maL1 i d r)+{-# INLINABLE alphaL1 #-}++autocorrL1 :: (Floating a, RealFloat a) => a -> a -> a -> a -> Fold a a+autocorrL1 i d maR corrR =+ case maL1 i d maR of+ (Fold maStep maBegin maDone) ->+ case corrL1 i d corrR of+ (Fold corrStep corrBegin corrDone) ->+ let begin = (maBegin, corrBegin)+ step (maAcc, corrAcc) a =+ ( maStep maAcc a+ , if isNaN (maDone maAcc)+ then corrAcc+ else corrStep corrAcc (maDone maAcc, a))+ done = corrDone . snd+ in Fold step begin done+{-# INLINABLE autocorrL1 #-}
src/Online/Quantiles.hs view
@@ -1,27 +1,32 @@ {-# LANGUAGE DataKinds #-}-{-# LANGUAGE ExistentialQuantification #-}-{-# LANGUAGE RankNTypes #-}-{-# LANGUAGE OverloadedStrings #-}-{-# LANGUAGE TypeSynonymInstances #-}-{-# OPTIONS_GHC -fno-warn-type-defaults #-} -module Online.Quantiles where+module Online.Quantiles+ ( tDigest+ , tDigestQuantiles+ , tDigestHist+ , OnlineTDigest(..)+ , onlineQuantiles+ , Online.Quantiles.median+ , onlineDigitize+ , onlineDigestHist+ )+where -import NumHask.Prelude import qualified Control.Foldl as L+import Data.List.NonEmpty (NonEmpty) import Data.TDigest import Data.TDigest.Internal-import qualified Data.Vector.Unboxed as VU+import Data.TDigest.Postprocess () import qualified Data.Vector.Algorithms.Heap as VHeap-import Data.List.NonEmpty (NonEmpty)-import Data.TDigest.Postprocess()+import qualified Data.Vector.Unboxed as VU+import NumHask.Prelude -- | a raw non-online tdigest fold tDigest :: L.Fold Double (TDigest 25) tDigest = L.Fold step begin done where step x a = insert a x- begin = tdigest ([]::[Double]) :: TDigest 25+ begin = tdigest ([] :: [Double]) :: TDigest 25 done = identity -- | non-online version@@ -29,7 +34,7 @@ tDigestQuantiles qs = L.Fold step begin done where step x a = insert a x- begin = tdigest ([]::[Double]) :: TDigest 25+ begin = tdigest ([] :: [Double]) :: TDigest 25 done x = fromMaybe nan . (`quantile` compress x) <$> qs -- | non-online version@@ -37,13 +42,17 @@ tDigestHist = L.Fold step begin done where step x a = insert a x- begin = tdigest ([]::[Double]) :: TDigest 25- done x = histogram . compress $ x+ begin = tdigest ([] :: [Double]) :: TDigest 25+ done = histogram . compress -data OnlineTDigest = OnlineTDigest { td :: TDigest 25, tdN :: Int, tdRate :: Double } deriving (Show)+data OnlineTDigest = OnlineTDigest+ { td :: TDigest 25+ , tdN :: Int+ , tdRate :: Double+ } deriving (Show) emptyOnlineTDigest :: Double -> OnlineTDigest-emptyOnlineTDigest r = OnlineTDigest (emptyTDigest :: TDigest n) 0 r+emptyOnlineTDigest = OnlineTDigest (emptyTDigest :: TDigest n) 0 -- | decaying quantiles based on the tdigest library onlineQuantiles :: Double -> [Double] -> L.Fold Double [Double]@@ -61,34 +70,39 @@ step x a = onlineInsert a x begin = emptyOnlineTDigest r done x = fromMaybe nan (quantile 0.5 t)- where- (OnlineTDigest t _ _) = onlineForceCompress x+ where+ (OnlineTDigest t _ _) = onlineForceCompress x onlineInsert' :: Double -> OnlineTDigest -> OnlineTDigest-onlineInsert' x (OnlineTDigest td n r) = OnlineTDigest (insertCentroid (x, r^^(-(fromIntegral $ n+1))) td) (n+1) r+onlineInsert' x (OnlineTDigest td' n r) =+ OnlineTDigest+ (insertCentroid (x, r ^^ (-(fromIntegral $ n + 1))) td')+ (n + 1)+ r onlineInsert :: Double -> OnlineTDigest -> OnlineTDigest onlineInsert x otd = onlineCompress (onlineInsert' x otd) onlineCompress :: OnlineTDigest -> OnlineTDigest-onlineCompress otd@(OnlineTDigest Nil _ _ ) = otd+onlineCompress otd@(OnlineTDigest Nil _ _) = otd onlineCompress otd@(OnlineTDigest t _ _)- | Data.TDigest.Internal.size t > relMaxSize * compression && Data.TDigest.Internal.size t > absMaxSize- = onlineForceCompress otd- | otherwise- = otd+ | Data.TDigest.Internal.size t > relMaxSize * compression &&+ Data.TDigest.Internal.size t > absMaxSize = onlineForceCompress otd+ | otherwise = otd where compression = 25 onlineForceCompress :: OnlineTDigest -> OnlineTDigest-onlineForceCompress otd@(OnlineTDigest Nil _ _ ) = otd+onlineForceCompress otd@(OnlineTDigest Nil _ _) = otd onlineForceCompress (OnlineTDigest t n r) = OnlineTDigest t' 0 r where- t' = NumHask.Prelude.foldl' (flip insertCentroid) emptyTDigest $- (\(m,w) -> (m, w*(r^^fromIntegral n))) . fst <$> VU.toList centroids+ t' =+ NumHask.Prelude.foldl' (flip insertCentroid) emptyTDigest $+ (\(m, w) -> (m, w * (r ^^ fromIntegral n))) . fst <$> VU.toList centroids -- Centroids are shuffled based on space centroids :: VU.Vector (Centroid, Double)- centroids = runST $ do+ centroids =+ runST $ do v <- toMVector t -- sort by cumulative weight VHeap.sortBy (comparing snd) v@@ -98,13 +112,19 @@ onlineDigitize :: Double -> [Double] -> L.Fold Double Int onlineDigitize r qs = L.Fold step begin done where- step (x,_) a = (onlineInsert a x, a)+ step (x, _) a = (onlineInsert a x, a) begin = (emptyOnlineTDigest r, nan)- done (x,l) = bucket' qs' l+ done (x, l) = bucket' qs' l where qs' = fromMaybe nan . (`quantile` t) <$> qs (OnlineTDigest t _ _) = onlineForceCompress x- bucket' xs l' = L.fold L.sum $ (\x' -> if x'>l' then 0 else 1) <$> xs+ bucket' xs l' =+ L.fold L.sum $+ (\x' ->+ if x' > l'+ then 0+ else 1) <$>+ xs -- | decaying histogram based on the tdigest library onlineDigestHist :: Double -> L.Fold Double (Maybe (NonEmpty HistBin))
− src/Online/Stats.hs
@@ -1,140 +0,0 @@-{-# LANGUAGE ExistentialQuantification #-}-{-# LANGUAGE GADTs #-}--module Online.Stats (- -- * convert a statistic to online- Averager,- online,-- -- * common statistics- av,-- -- * online statistics- ma,- absma,- sqma,- std,- cov,- corr,- corrGauss,- beta,- alpha,- autocorr- ) where--import Protolude-import qualified Control.Foldl as L-import Control.Foldl (Fold(..))---- | Most common statistics are averages.-newtype Averager a b = Averager { _averager :: (a, b)}--instance (Monoid a, Monoid b) => Monoid (Averager a b) where- mempty = Averager (mempty, mempty)- mappend (Averager (s,c)) (Averager (s',c')) = Averager (mappend s s', mappend c c')---- | online takes a function and turns it into a `Control.Foldl.Fold` where the step is an incremental update of the (isomorphic) statistic.-online :: (Fractional b) => (a -> b) -> (b -> b) -> Fold a b-online f g = Fold step begin extract- where- begin = Averager (0, 0)- step (Averager (s,c)) a = Averager (g $ s+f a,g $ c+1)- extract (Averager (s,c)) = s/c---- | average-av :: (Fractional a) => Fold a a-av = Fold step begin extract- where- begin = Averager (0, 0)- step (Averager (s,c)) a = Averager (s+a,c+1)- extract (Averager (s,c)) = s/c-{-# INLINABLE av #-}---- | moving average-ma :: (Fractional a) => a -> Fold a a-ma r = online identity (* r)-{-# INLINABLE ma #-}---- | absolute average-absma :: (Fractional a) => a -> Fold a a-absma r = online abs (* r)-{-# INLINABLE absma #-}---- | average square-sqma :: (Fractional a) => a -> Fold a a-sqma r = online (\x -> x*x) (* r)-{-# INLINABLE sqma #-}---- | standard deviation-std :: (Floating a) => a -> Fold a a-std r = (\s ss -> sqrt (ss - s**2)) <$> ma r <*> sqma r-{-# INLINABLE std #-}---- | the covariance of a tuple--- given an underlying central tendency fold-cov :: (Floating a) => Fold a a -> Fold (a,a) a-cov m =- (\xy x' y' -> xy - x' * y') <$>- L.premap (uncurry (*)) m <*>- L.premap fst m <*>- L.premap snd m-{-# INLINABLE cov #-}---- | correlation of a tuple, specialised to Guassian-corrGauss :: (Floating a) => a -> Fold (a,a) a-corrGauss r = (\cov' stdx stdy -> cov' / (stdx * stdy)) <$> cov (ma r) <*> L.premap fst (std r) <*> L.premap snd (std r)-{-# INLINABLE corrGauss #-}---- | a generalised version of correlation of a tuple-corr :: (Floating a) => Fold a a -> Fold a a -> Fold (a,a) a-corr central deviation =- (\cov' stdx stdy -> cov' / (stdx * stdy)) <$>- cov central <*>- L.premap fst deviation <*>- L.premap snd deviation-{-# INLINABLE corr #-}---- | the beta in a simple linear regression of a tuple--- given an underlying central tendency fold-beta :: (Floating a) => Fold a a -> Fold (a,a) a-beta m =- (\xy x' y' x2 -> (xy - x'*y')/(x2 - x'*x')) <$>- L.premap (uncurry (*)) m <*>- L.premap fst m <*>- L.premap snd m <*>- L.premap (\(x,_) -> x*x) m-{-# INLINABLE beta #-}---- | the alpha of a tuple-alpha :: (Floating a) => Fold a a -> Fold (a,a) a-alpha m =- (\y b x -> y - b * x) <$>- L.premap fst m <*> beta m <*> L.premap snd m-{-# INLINABLE alpha #-}--{-| autocorrelation is a slippery concept. This method starts with the concept that there is an underlying random error process (e), and autocorrelation is a process on top of that ie for a one-step correlation relationship.--value@t = e@t + k * e@t-1--where k is the autocorrelation.--There are thus two online rates needed: one for the average being considered to be the dependent variable, and one for the online of the correlation calculation between the most recent value and the moving average. For example,-->>> L.fold (autocorr 0 1)--would estimate the one-step autocorrelation relationship of the previous value and the current value over the entire sample set. ---}-autocorr :: (Floating a, RealFloat a) => Fold a a -> Fold (a,a) a -> Fold a a-autocorr central corrf =- case central of- (Fold mStep mBegin mDone) ->- case corrf of- (Fold dStep dBegin dDone) ->- let begin = (mBegin, dBegin)- step (mAcc,dAcc) a = (mStep mAcc a,- if isNaN (mDone mAcc)- then dAcc- else dStep dAcc (mDone mAcc, a))- done = dDone . snd in- Fold step begin done
− src/Online/StatsL1.hs
@@ -1,112 +0,0 @@-{-# LANGUAGE ExistentialQuantification #-}-{-# LANGUAGE GADTs #-}--module Online.StatsL1 (- -- * convert a statistic to an online median stat equivalent to L1- Medianer(..),- onlineL1,- onlineL1',-- -- * common statistics- avL1,-- -- * online statistics- maL1,- absmaL1,- covL1,- corrL1,- betaL1,- alphaL1,- autocorrL1- ) where--import Protolude-import qualified Control.Foldl as L-import Control.Foldl (Fold(..))---- | A rough Median.--- The average absolute value of the stat is used to callibrate estimate drift towards the medium-data Medianer a b = Medianer { medAbsSum :: a, medCount :: b, medianEst :: a}---- | onlineL1' takes a function and turns it into a `Control.Foldl.Fold` where the step is an incremental update of an (isomorphic) median statistic.-onlineL1' :: (Ord b, Fractional b) => b -> b -> (a -> b) -> (b -> b) -> Fold a (b,b)-onlineL1' i d f g = Fold step begin extract- where- begin = Medianer 0 0 0- step (Medianer s c m) a = Medianer (g $ s+abs (f a)) (g $ c+1) ((1-d) * (m + s' * i * s/c') + d * f a)- where- c' = if c == 0 then 1 else c- s'- | f a > m = 1- | f a < m = -1- | otherwise = 0- extract (Medianer s c m) = (s/c,m)---- | onlineL1 takes a function and turns it into a `Control.Foldl.Fold` where the step is an incremental update of an (isomorphic) median statistic.-onlineL1 :: (Ord b, Fractional b) => b -> b -> (a -> b) -> (b -> b) -> Fold a b-onlineL1 i d f g = snd <$> onlineL1' i d f g----- | averageL1-avL1 :: (Ord a, Fractional a) => a -> Fold a a-avL1 i = Fold step begin extract- where- begin = Medianer 0 0 0- step (Medianer s c m) a = Medianer (s+a) (c+1) (m + s' * i * s/c)- where- s'- | a > m = 1- | a < m = -1- | otherwise = 0- extract (Medianer _ _ m) = m-{-# INLINABLE avL1 #-}---- | moving median-maL1 :: (Ord a, Fractional a) => a -> a -> a -> Fold a a-maL1 i d r = onlineL1 i d identity (*r)-{-# INLINABLE maL1 #-}---- | moving absolute deviation-absmaL1 :: (Ord a, Fractional a) => a -> a -> a -> Fold a a-absmaL1 i d r = fst <$> onlineL1' i d identity (*r)-{-# INLINABLE absmaL1 #-}---- | covariance of a tuple-covL1 :: (Ord a, Fractional a) => a -> a -> a -> Fold (a,a) a-covL1 i d r = (\xy xbar ybar -> xy - xbar * ybar) <$> onlineL1 i d (uncurry (*)) (*r) <*> onlineL1 i d fst (*r) <*> onlineL1 i d snd (*r)-{-# INLINABLE covL1 #-}---- | correlation of a tuple-corrL1 :: (Ord a, Floating a) => a -> a -> a -> Fold (a,a) a-corrL1 i d r = (\cov' stdx stdy -> cov' / (stdx * stdy)) <$> covL1 i d r <*> L.premap fst (absmaL1 i d r) <*> L.premap snd (absmaL1 i d r)-{-# INLINABLE corrL1 #-}---- | the beta in a simple linear regression of a tuple-betaL1 :: (Ord a, Floating a) => a -> a -> a -> Fold (a,a) a-betaL1 i d r =- (\xy x' y' x2 -> (xy - x'*y')/(x2 - x'*x')) <$>- L.premap (uncurry (*)) (maL1 i d r) <*>- L.premap fst (maL1 i d r) <*>- L.premap snd (maL1 i d r) <*>- L.premap (\(x,_) -> x*x) (maL1 i d r)-{-# INLINABLE betaL1 #-}---- | the alpha in a simple linear regression of `snd` on `fst`-alphaL1 :: (Ord a, Floating a) => a -> a -> a -> Fold (a,a) a-alphaL1 i d r = (\y b x -> y - b * x) <$> L.premap fst (maL1 i d r) <*> betaL1 i d r <*> L.premap snd (maL1 i d r)-{-# INLINABLE alphaL1 #-}--autocorrL1 :: (Floating a, RealFloat a) => a -> a -> a -> a -> Fold a a-autocorrL1 i d maR corrR = - case maL1 i d maR of- (Fold maStep maBegin maDone) ->- case corrL1 i d corrR of- (Fold corrStep corrBegin corrDone) ->- let begin = (maBegin, corrBegin)- step (maAcc,corrAcc) a = (maStep maAcc a,- if isNaN (maDone maAcc)- then corrAcc- else corrStep corrAcc (maDone maAcc, a)) - done = corrDone . snd in- Fold step begin done-{-# INLINABLE autocorrL1 #-}
stack.yaml view
@@ -1,7 +1,6 @@-resolver: lts-8.23+resolver: nightly-2017-11-22 packages: - '.' -extra-deps:- - numhask-0.0.7+extra-deps: []
+ test/test.hs view
@@ -0,0 +1,19 @@+{-# OPTIONS_GHC -Wall #-}++module Main where++import Protolude+import Test.Tasty (TestTree, testGroup, defaultMain)+import Test.DocTest++main :: IO ()+main = do+ doctest ["src/Online/Averages.hs", "src/Online/Medians.hs"]+ defaultMain tests++tests :: TestTree+tests =+ testGroup ""+ [+ ]+