{-# OPTIONS_GHC -fplugin=GHC.TypeLits.KnownNat.Solver -fplugin=GHC.TypeLits.Normalise -fconstraint-solver-iterations=10 #-}
{-# LANGUAGE
RankNTypes,
TypeOperators,
FlexibleContexts,
ScopedTypeVariables
#-}
module Goal.Probability
(
module Goal.Probability.Statistical
, module Goal.Probability.ExponentialFamily
, module Goal.Probability.Conditional
, module Goal.Probability.Distributions
, module Goal.Probability.Distributions.Gaussian
, module Goal.Probability.Distributions.CoMPoisson
, shuffleList
, resampleVector
, subsampleVector
, noisyFunction
, minibatcher
, estimateMeanVariance
, estimateMeanSD
, estimateFanoFactor
, estimateCoefficientOfVariation
, estimateCorrelation
, estimateCorrelations
, histograms
, akaikesInformationCriterion
, bayesianInformationCriterion
) where
import Goal.Probability.Statistical
import Goal.Probability.ExponentialFamily
import Goal.Probability.Conditional
import Goal.Probability.Distributions
import Goal.Probability.Distributions.Gaussian
import Goal.Probability.Distributions.CoMPoisson
import Goal.Core
import Goal.Geometry
import qualified Goal.Core.Vector.Boxed as B
import qualified Goal.Core.Vector.Storable as S
import qualified Goal.Core.Vector.Generic.Mutable as M
import qualified Goal.Core.Vector.Generic as G
import qualified Data.Vector.Generic.Mutable.Base as MV
import qualified Data.Vector as V
import qualified Statistics.Sample as STAT hiding (range)
import qualified Statistics.Sample.Histogram as STAT
import qualified Data.Vector.Storable as VS
import qualified System.Random.MWC as R
import qualified System.Random.MWC.Distributions as R
estimateMeanVariance
:: Traversable f
=> f Double
-> (Double,Double)
estimateMeanVariance :: f Double -> (Double, Double)
estimateMeanVariance f Double
xs = Vector Double -> (Double, Double)
forall (v :: Type -> Type).
Vector v Double =>
v Double -> (Double, Double)
STAT.meanVarianceUnb (Vector Double -> (Double, Double))
-> ([Double] -> Vector Double) -> [Double] -> (Double, Double)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
VS.fromList ([Double] -> (Double, Double)) -> [Double] -> (Double, Double)
forall a b. (a -> b) -> a -> b
$ f Double -> [Double]
forall (t :: Type -> Type) a. Foldable t => t a -> [a]
toList f Double
xs
estimateMeanSD
:: Traversable f
=> f Double
-> (Double,Double)
estimateMeanSD :: f Double -> (Double, Double)
estimateMeanSD f Double
xs =
let (Double
mu,Double
vr) = f Double -> (Double, Double)
forall (f :: Type -> Type).
Traversable f =>
f Double -> (Double, Double)
estimateMeanVariance f Double
xs
in (Double
mu,Double -> Double
forall a. Floating a => a -> a
sqrt Double
vr)
estimateFanoFactor
:: Traversable f
=> f Double
-> Double
estimateFanoFactor :: f Double -> Double
estimateFanoFactor f Double
xs =
let (Double
mu,Double
vr) = f Double -> (Double, Double)
forall (f :: Type -> Type).
Traversable f =>
f Double -> (Double, Double)
estimateMeanVariance f Double
xs
in Double
vr Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
mu
estimateCoefficientOfVariation :: Traversable f => f Double -> Double
estimateCoefficientOfVariation :: f Double -> Double
estimateCoefficientOfVariation f Double
zs =
let (Double
mu,Double
vr) = f Double -> (Double, Double)
forall (f :: Type -> Type).
Traversable f =>
f Double -> (Double, Double)
estimateMeanVariance f Double
zs
in Double -> Double
forall a. Floating a => a -> a
sqrt Double
vr Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
mu
estimateCorrelations
:: forall k x v . (G.VectorClass v x, G.VectorClass v Double, KnownNat k, Real x)
=> [G.Vector v k x]
-> S.Matrix k k Double
estimateCorrelations :: [Vector v k x] -> Matrix k k Double
estimateCorrelations [Vector v k x]
zs =
let mnrm :: Source # MultivariateNormal k
mnrm :: Source # MultivariateNormal k
mnrm = Sample (MultivariateNormal k) -> Source # MultivariateNormal k
forall c x. MaximumLikelihood c x => Sample x -> c # x
mle (Sample (MultivariateNormal k) -> Source # MultivariateNormal k)
-> Sample (MultivariateNormal k) -> Source # MultivariateNormal k
forall a b. (a -> b) -> a -> b
$ Vector v k Double -> Vector Vector k Double
forall (v :: Type -> Type) a (w :: Type -> Type) (n :: Nat).
(Vector v a, Vector w a) =>
Vector v n a -> Vector w n a
G.convert (Vector v k Double -> Vector Vector k Double)
-> (Vector v k x -> Vector v k Double)
-> Vector v k x
-> Vector Vector k Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (x -> Double) -> Vector v k x -> Vector v k Double
forall (v :: Type -> Type) a b (n :: Nat).
(Vector v a, Vector v b) =>
(a -> b) -> Vector v n a -> Vector v n b
G.map x -> Double
forall a b. (Real a, Fractional b) => a -> b
realToFrac (Vector v k x -> Vector Vector k Double)
-> [Vector v k x] -> [Vector Vector k Double]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector v k x]
zs
in (Source # MultivariateNormal k) -> Matrix k k Double
forall (k :: Nat).
KnownNat k =>
(Source # MultivariateNormal k) -> Matrix k k Double
multivariateNormalCorrelations Source # MultivariateNormal k
mnrm
estimateCorrelation
:: [(Double,Double)]
-> Double
estimateCorrelation :: [(Double, Double)] -> Double
estimateCorrelation [(Double, Double)]
zs = Vector (Double, Double) -> Double
forall (v :: Type -> Type).
(Vector v (Double, Double), Vector v Double) =>
v (Double, Double) -> Double
STAT.correlation (Vector (Double, Double) -> Double)
-> Vector (Double, Double) -> Double
forall a b. (a -> b) -> a -> b
$ [(Double, Double)] -> Vector (Double, Double)
forall a. [a] -> Vector a
V.fromList [(Double, Double)]
zs
histograms
:: Int
-> Maybe (Double, Double)
-> [[Double]]
-> ([Double],[[Int]],[[Double]])
histograms :: Int
-> Maybe (Double, Double)
-> [[Double]]
-> ([Double], [[Int]], [[Double]])
histograms Int
nbns Maybe (Double, Double)
mmnmx [[Double]]
smpss =
let (Double
mn,Double
mx) = case Maybe (Double, Double)
mmnmx of
Just (Double
mn0,Double
mx0) -> (Double
mn0,Double
mx0)
Maybe (Double, Double)
Nothing -> Int -> Vector Double -> (Double, Double)
forall (v :: Type -> Type).
Vector v Double =>
Int -> v Double -> (Double, Double)
STAT.range Int
nbns (Vector Double -> (Double, Double))
-> ([Double] -> Vector Double) -> [Double] -> (Double, Double)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
VS.fromList ([Double] -> (Double, Double)) -> [Double] -> (Double, Double)
forall a b. (a -> b) -> a -> b
$ [[Double]] -> [Double]
forall (t :: Type -> Type) a. Foldable t => t [a] -> [a]
concat [[Double]]
smpss
stp :: Double
stp = (Double
mx Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
mn) Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
nbns
bns :: [Double]
bns = Int -> [Double] -> [Double]
forall a. Int -> [a] -> [a]
take Int
nbns [ Double
mn Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
stpDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
stp Double -> Double -> Double
forall a. Num a => a -> a -> a
* Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n | Int
n <- [Int
0 :: Int,Int
1..] ]
hsts :: [[Int]]
hsts = Vector Int -> [Int]
forall a. Storable a => Vector a -> [a]
VS.toList (Vector Int -> [Int])
-> ([Double] -> Vector Int) -> [Double] -> [Int]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Double -> Double -> Vector Double -> Vector Int
forall b a (v0 :: Type -> Type) (v1 :: Type -> Type).
(Num b, RealFrac a, Vector v0 a, Vector v1 b) =>
Int -> a -> a -> v0 a -> v1 b
STAT.histogram_ Int
nbns Double
mn Double
mx (Vector Double -> Vector Int)
-> ([Double] -> Vector Double) -> [Double] -> Vector Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
VS.fromList ([Double] -> [Int]) -> [[Double]] -> [[Int]]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [[Double]]
smpss
ttls :: [Int]
ttls = [Int] -> Int
forall (t :: Type -> Type) a. (Foldable t, Num a) => t a -> a
sum ([Int] -> Int) -> [[Int]] -> [Int]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [[Int]]
hsts
dnss :: [[Double]]
dnss = do
([Int]
hst,Int
ttl) <- [[Int]] -> [Int] -> [([Int], Int)]
forall a b. [a] -> [b] -> [(a, b)]
zip [[Int]]
hsts [Int]
ttls
[Double] -> [[Double]]
forall (m :: Type -> Type) a. Monad m => a -> m a
return ([Double] -> [[Double]]) -> [Double] -> [[Double]]
forall a b. (a -> b) -> a -> b
$ if Int
ttl Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0
then []
else (Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/(Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ttl Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
stp)) (Double -> Double) -> (Int -> Double) -> Int -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Double) -> [Int] -> [Double]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Int]
hst
in ([Double]
bns,[[Int]]
hsts,[[Double]]
dnss)
shuffleList :: [a] -> Random [a]
shuffleList :: [a] -> Random [a]
shuffleList [a]
xs = Vector a -> [a]
forall a. Vector a -> [a]
V.toList (Vector a -> [a]) -> Random (Vector a) -> Random [a]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> (forall s. Gen s -> ST s (Vector a)) -> Random (Vector a)
forall a. (forall s. Gen s -> ST s a) -> Random a
Random (Vector a -> Gen s -> ST s (Vector a)
forall g (m :: Type -> Type) (v :: Type -> Type) a.
(StatefulGen g m, PrimMonad m, Vector v a) =>
v a -> g -> m (v a)
R.uniformShuffle ([a] -> Vector a
forall a. [a] -> Vector a
V.fromList [a]
xs))
minibatcher :: Int -> [x] -> Chain Random [x]
minibatcher :: Int -> [x] -> Chain Random [x]
minibatcher Int
nbtch [x]
xs0 = [x] -> (() -> [x] -> Random ([x], [x])) -> Chain Random [x]
forall (m :: Type -> Type) acc a b.
Monad m =>
acc -> (a -> acc -> m (b, acc)) -> Circuit m a b
accumulateFunction [] ((() -> [x] -> Random ([x], [x])) -> Chain Random [x])
-> (() -> [x] -> Random ([x], [x])) -> Chain Random [x]
forall a b. (a -> b) -> a -> b
$ \() [x]
xs ->
if [x] -> Int
forall (t :: Type -> Type) a. Foldable t => t a -> Int
length (Int -> [x] -> [x]
forall a. Int -> [a] -> [a]
take Int
nbtch [x]
xs) Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
nbtch
then do
[x]
xs1 <- [x] -> Random [x]
forall a. [a] -> Random [a]
shuffleList [x]
xs0
let ([x]
hds',[x]
tls') = Int -> [x] -> ([x], [x])
forall a. Int -> [a] -> ([a], [a])
splitAt Int
nbtch ([x]
xs [x] -> [x] -> [x]
forall a. [a] -> [a] -> [a]
++ [x]
xs1)
([x], [x]) -> Random ([x], [x])
forall (m :: Type -> Type) a. Monad m => a -> m a
return ([x]
hds',[x]
tls')
else do
let ([x]
hds',[x]
tls') = Int -> [x] -> ([x], [x])
forall a. Int -> [a] -> ([a], [a])
splitAt Int
nbtch [x]
xs
([x], [x]) -> Random ([x], [x])
forall (m :: Type -> Type) a. Monad m => a -> m a
return ([x]
hds',[x]
tls')
resampleVector :: (KnownNat n, KnownNat k) => B.Vector n x -> Random (B.Vector k x)
resampleVector :: Vector n x -> Random (Vector k x)
resampleVector Vector n x
xs = do
Vector k Int
ks <- Random Int -> Random (Vector k Int)
forall (n :: Nat) (m :: Type -> Type) a.
(KnownNat n, Monad m) =>
m a -> m (Vector n a)
B.replicateM (Random Int -> Random (Vector k Int))
-> Random Int -> Random (Vector k Int)
forall a b. (a -> b) -> a -> b
$ (forall s. Gen s -> ST s Int) -> Random Int
forall a. (forall s. Gen s -> ST s a) -> Random a
Random ((Int, Int) -> Gen (PrimState (ST s)) -> ST s Int
forall a (m :: Type -> Type).
(Variate a, PrimMonad m) =>
(a, a) -> Gen (PrimState m) -> m a
R.uniformR (Int
0, Vector n x -> Int
forall (n :: Nat) a. KnownNat n => Vector n a -> Int
B.length Vector n x
xsInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1))
Vector k x -> Random (Vector k x)
forall (m :: Type -> Type) a. Monad m => a -> m a
return (Vector k x -> Random (Vector k x))
-> Vector k x -> Random (Vector k x)
forall a b. (a -> b) -> a -> b
$ Vector n x -> Vector k Int -> Vector k x
forall (m :: Nat) a (n :: Nat).
Vector m a -> Vector n Int -> Vector n a
B.backpermute Vector n x
xs Vector k Int
ks
noisyFunction
:: (Generative c x, Num (SamplePoint x))
=> Point c x
-> (y -> SamplePoint x)
-> y
-> Random (SamplePoint x)
noisyFunction :: Point c x -> (y -> SamplePoint x) -> y -> Random (SamplePoint x)
noisyFunction Point c x
m y -> SamplePoint x
f y
x = do
SamplePoint x
ns <- Point c x -> Random (SamplePoint x)
forall c x. Generative c x => Point c x -> Random (SamplePoint x)
samplePoint Point c x
m
SamplePoint x -> Random (SamplePoint x)
forall (m :: Type -> Type) a. Monad m => a -> m a
return (SamplePoint x -> Random (SamplePoint x))
-> SamplePoint x -> Random (SamplePoint x)
forall a b. (a -> b) -> a -> b
$ y -> SamplePoint x
f y
x SamplePoint x -> SamplePoint x -> SamplePoint x
forall a. Num a => a -> a -> a
+ SamplePoint x
ns
subsampleVector
:: forall k m v x . (KnownNat k, KnownNat m, G.VectorClass v x)
=> G.Vector v (k + m) x
-> Random (G.Vector v k x)
subsampleVector :: Vector v (k + m) x -> Random (Vector v k x)
subsampleVector Vector v (k + m) x
v = (forall s. Gen s -> ST s (Vector v k x)) -> Random (Vector v k x)
forall a. (forall s. Gen s -> ST s a) -> Random a
Random ((forall s. Gen s -> ST s (Vector v k x)) -> Random (Vector v k x))
-> (forall s. Gen s -> ST s (Vector v k x))
-> Random (Vector v k x)
forall a b. (a -> b) -> a -> b
$ \Gen s
gn -> do
let k :: Int
k = Proxy k -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy k
forall k (t :: k). Proxy t
Proxy :: Proxy k)
MVector (Mutable v) (k + m) s x
mv <- Vector v (k + m) x
-> ST s (MVector (Mutable v) (k + m) (PrimState (ST s)) x)
forall (m :: Type -> Type) (v :: Type -> Type) a (n :: Nat).
(PrimMonad m, Vector v a) =>
Vector v n a -> m (MVector (Mutable v) n (PrimState m) a)
G.thaw Vector v (k + m) x
v
Int
-> MVector (Mutable v) (k + m) (PrimState (ST s)) x
-> Gen (PrimState (ST s))
-> ST s ()
forall (n :: Nat) (m :: Type -> Type) (v :: Type -> Type -> Type)
a.
(KnownNat n, PrimMonad m, MVector v a) =>
Int -> MVector v n (PrimState m) a -> Gen (PrimState m) -> m ()
randomSubSample0 Int
k MVector (Mutable v) (k + m) s x
MVector (Mutable v) (k + m) (PrimState (ST s)) x
mv Gen s
Gen (PrimState (ST s))
gn
Vector v (k + m) x
v' <- MVector (Mutable v) (k + m) (PrimState (ST s)) x
-> ST s (Vector v (k + m) x)
forall (m :: Type -> Type) (v :: Type -> Type) a (n :: Nat).
(PrimMonad m, Vector v a) =>
MVector (Mutable v) n (PrimState m) a -> m (Vector v n a)
G.unsafeFreeze MVector (Mutable v) (k + m) s x
MVector (Mutable v) (k + m) (PrimState (ST s)) x
mv
let foo :: (G.Vector v k x, G.Vector v m x)
foo :: (Vector v k x, Vector v m x)
foo = Vector v (k + m) x -> (Vector v k x, Vector v m x)
forall (v :: Type -> Type) (n :: Nat) (m :: Nat) a.
(KnownNat n, Vector v a) =>
Vector v (n + m) a -> (Vector v n a, Vector v m a)
G.splitAt Vector v (k + m) x
v'
Vector v k x -> ST s (Vector v k x)
forall (m :: Type -> Type) a. Monad m => a -> m a
return (Vector v k x -> ST s (Vector v k x))
-> Vector v k x -> ST s (Vector v k x)
forall a b. (a -> b) -> a -> b
$ (Vector v k x, Vector v m x) -> Vector v k x
forall a b. (a, b) -> a
fst (Vector v k x, Vector v m x)
foo
randomSubSample0
:: (KnownNat n, PrimMonad m, MV.MVector v a)
=> Int -> G.MVector v n (PrimState m) a -> R.Gen (PrimState m) -> m ()
randomSubSample0 :: Int -> MVector v n (PrimState m) a -> Gen (PrimState m) -> m ()
randomSubSample0 Int
k MVector v n (PrimState m) a
v Gen (PrimState m)
gn = Int -> m ()
looper Int
0
where n :: Int
n = MVector v n (PrimState m) a -> Int
forall (v :: Type -> Type -> Type) (n :: Nat) s a.
KnownNat n =>
MVector v n s a -> Int
M.length MVector v n (PrimState m) a
v
looper :: Int -> m ()
looper Int
i
| Int
i Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
k = () -> m ()
forall (m :: Type -> Type) a. Monad m => a -> m a
return ()
| Bool
otherwise = do
Int
j <- (Int, Int) -> Gen (PrimState m) -> m Int
forall a (m :: Type -> Type).
(Variate a, PrimMonad m) =>
(a, a) -> Gen (PrimState m) -> m a
R.uniformR (Int
i,Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1) Gen (PrimState m)
gn
MVector v n (PrimState m) a -> Int -> Int -> m ()
forall (v :: Type -> Type -> Type) (n :: Nat) (m :: Type -> Type)
a.
(PrimMonad m, MVector v a) =>
MVector v n (PrimState m) a -> Int -> Int -> m ()
M.unsafeSwap MVector v n (PrimState m) a
v Int
i Int
j
Int -> m ()
looper (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
akaikesInformationCriterion
:: forall c x s . (Manifold x, LogLikelihood c x s)
=> c # x
-> [s]
-> Double
akaikesInformationCriterion :: (c # x) -> [s] -> Double
akaikesInformationCriterion c # x
p [s]
xs =
let d :: Natural
d = Proxy (Dimension x) -> Natural
forall (n :: Nat) (proxy :: Nat -> Type).
KnownNat n =>
proxy n -> Natural
natVal (Proxy (Dimension x)
forall k (t :: k). Proxy t
Proxy :: Proxy (Dimension x))
in Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Natural -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Natural
d Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* [s] -> (c # x) -> Double
forall c x s. LogLikelihood c x s => [s] -> (c # x) -> Double
logLikelihood [s]
xs c # x
p Double -> Double -> Double
forall a. Num a => a -> a -> a
* Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral ([s] -> Int
forall (t :: Type -> Type) a. Foldable t => t a -> Int
length [s]
xs)
bayesianInformationCriterion
:: forall c x s . (LogLikelihood c x s, Manifold x)
=> c # x
-> [s]
-> Double
bayesianInformationCriterion :: (c # x) -> [s] -> Double
bayesianInformationCriterion c # x
p [s]
xs =
let d :: Natural
d = Proxy (Dimension x) -> Natural
forall (n :: Nat) (proxy :: Nat -> Type).
KnownNat n =>
proxy n -> Natural
natVal (Proxy (Dimension x)
forall k (t :: k). Proxy t
Proxy :: Proxy (Dimension x))
n :: Int
n = [s] -> Int
forall (t :: Type -> Type) a. Foldable t => t a -> Int
length [s]
xs
in Double -> Double
forall a. Floating a => a -> a
log (Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n) Double -> Double -> Double
forall a. Num a => a -> a -> a
* Natural -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Natural
d Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* [s] -> (c # x) -> Double
forall c x s. LogLikelihood c x s => [s] -> (c # x) -> Double
logLikelihood [s]
xs c # x
p Double -> Double -> Double
forall a. Num a => a -> a -> a
* Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral ([s] -> Int
forall (t :: Type -> Type) a. Foldable t => t a -> Int
length [s]
xs)