module MCSP.Algorithms.Vector (
Default,
zeros,
replicate,
map,
choose,
(.+),
(.-),
(.*),
(.*.),
sum,
sort,
sortOn,
sortBy,
argSort,
sortLike,
normalized,
standardized,
sumM,
replicateM,
uniformN,
uniformSN,
uniformRN,
weighted,
weightedN,
choice,
) where
import Control.Applicative (pure)
import Control.Exception.Extra (errorWithoutStackTrace)
import Control.Monad (Monad, join, sequence)
import Data.Bool (Bool (..), bool, otherwise, (||))
import Data.Eq (Eq (..))
import Data.Foldable1 (foldl1')
import Data.Function (id, on, ($), (.))
import Data.Functor ((<$>))
import Data.Int (Int)
import Data.List.NonEmpty (NonEmpty)
import Data.Ord (Ord (..), Ordering)
import Data.Vector.Algorithms.Merge qualified as Vector (sort, sortBy)
import Data.Vector.Generic qualified as Vector (maximumOn, sum)
import Data.Vector.Unboxed (
Unbox,
Vector,
create,
length,
map,
modify,
null,
replicate,
replicateM,
unsafeBackpermute,
unsafeIndex,
zipWith,
)
import Data.Vector.Unboxed.Mutable (generate)
import GHC.Float (Double, Float, Floating, sqrt)
import GHC.Num (Num (..))
import GHC.Real (Fractional, fromIntegral, (/))
import Text.Printf (printf)
import MCSP.System.Random (Random, Variate, uniformR, weightedChoice)
type Default = Float
withSameLength :: (Unbox a, Unbox b) => (Vector a -> Vector b -> c) -> Vector a -> Vector b -> c
withSameLength :: forall a b c.
(Unbox a, Unbox b) =>
(Vector a -> Vector b -> c) -> Vector a -> Vector b -> c
withSameLength Vector a -> Vector b -> c
f Vector a
v1 Vector b
v2 =
if Vector a -> Int
forall a. Unbox a => Vector a -> Int
length Vector a
v1 Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Vector b -> Int
forall a. Unbox a => Vector a -> Int
length Vector b
v2
then Vector a -> Vector b -> c
f Vector a
v1 Vector b
v2
else [Char] -> c
forall a. [Char] -> a
errorWithoutStackTrace ([Char] -> c) -> [Char] -> c
forall a b. (a -> b) -> a -> b
$ [Char] -> Int -> Int -> [Char]
forall r. PrintfType r => [Char] -> r
printf [Char]
"length mismatch: %d != %d" (Vector a -> Int
forall a. Unbox a => Vector a -> Int
length Vector a
v1) (Vector b -> Int
forall a. Unbox a => Vector a -> Int
length Vector b
v2)
{-# INLINE withSameLength #-}
zeros :: (Unbox a, Num a) => Int -> Vector a
zeros :: forall a. (Unbox a, Num a) => Int -> Vector a
zeros Int
s = Int -> a -> Vector a
forall a. Unbox a => Int -> a -> Vector a
replicate Int
s a
0
{-# SPECIALIZE zeros :: Int -> Vector Default #-}
choose :: Unbox a => a -> a -> Vector Bool -> Vector a
choose :: forall a. Unbox a => a -> a -> Vector Bool -> Vector a
choose a
falsy a
truthy = (Bool -> a) -> Vector Bool -> Vector a
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
map (a -> a -> Bool -> a
forall a. a -> a -> Bool -> a
bool a
falsy a
truthy)
{-# SPECIALIZE choose :: Default -> Default -> Vector Bool -> Vector Default #-}
infixl 7 .*, .*.
infixl 6 .+, .-
(.+) :: (Unbox a, Num a) => Vector a -> Vector a -> Vector a
.+ :: forall a. (Unbox a, Num a) => Vector a -> Vector a -> Vector a
(.+) = (a -> a -> a) -> Vector a -> Vector a -> Vector a
forall a b c.
(Unbox a, Unbox b, Unbox c) =>
(a -> b -> c) -> Vector a -> Vector b -> Vector c
zipWith a -> a -> a
forall a. Num a => a -> a -> a
(+)
{-# SPECIALIZE (.+) :: Vector Default -> Vector Default -> Vector Default #-}
(.-) :: (Unbox a, Num a) => Vector a -> Vector a -> Vector a
.- :: forall a. (Unbox a, Num a) => Vector a -> Vector a -> Vector a
(.-) = (Vector a -> Vector a -> Vector a)
-> Vector a -> Vector a -> Vector a
forall a b c.
(Unbox a, Unbox b) =>
(Vector a -> Vector b -> c) -> Vector a -> Vector b -> c
withSameLength ((Vector a -> Vector a -> Vector a)
-> Vector a -> Vector a -> Vector a)
-> (Vector a -> Vector a -> Vector a)
-> Vector a
-> Vector a
-> Vector a
forall a b. (a -> b) -> a -> b
$ (a -> a -> a) -> Vector a -> Vector a -> Vector a
forall a b c.
(Unbox a, Unbox b, Unbox c) =>
(a -> b -> c) -> Vector a -> Vector b -> Vector c
zipWith (-)
{-# SPECIALIZE (.-) :: Vector Default -> Vector Default -> Vector Default #-}
(.*) :: (Unbox a, Num a) => Vector a -> Vector a -> Vector a
.* :: forall a. (Unbox a, Num a) => Vector a -> Vector a -> Vector a
(.*) = (Vector a -> Vector a -> Vector a)
-> Vector a -> Vector a -> Vector a
forall a b c.
(Unbox a, Unbox b) =>
(Vector a -> Vector b -> c) -> Vector a -> Vector b -> c
withSameLength ((Vector a -> Vector a -> Vector a)
-> Vector a -> Vector a -> Vector a)
-> (Vector a -> Vector a -> Vector a)
-> Vector a
-> Vector a
-> Vector a
forall a b. (a -> b) -> a -> b
$ (a -> a -> a) -> Vector a -> Vector a -> Vector a
forall a b c.
(Unbox a, Unbox b, Unbox c) =>
(a -> b -> c) -> Vector a -> Vector b -> Vector c
zipWith a -> a -> a
forall a. Num a => a -> a -> a
(*)
{-# SPECIALIZE (.*) :: Vector Default -> Vector Default -> Vector Default #-}
(.*.) :: (Unbox a, Num a) => a -> Vector a -> Vector a
a
factor .*. :: forall a. (Unbox a, Num a) => a -> Vector a -> Vector a
.*. Vector a
vector = (a -> a) -> Vector a -> Vector a
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
map (a
factor *) Vector a
vector
{-# SPECIALIZE (.*.) :: Default -> Vector Default -> Vector Default #-}
sum :: (Unbox a, Num a) => NonEmpty (Vector a) -> Vector a
sum :: forall a. (Unbox a, Num a) => NonEmpty (Vector a) -> Vector a
sum = (Vector a -> Vector a -> Vector a)
-> NonEmpty (Vector a) -> Vector a
forall (t :: * -> *) a. Foldable1 t => (a -> a -> a) -> t a -> a
foldl1' Vector a -> Vector a -> Vector a
forall a. (Unbox a, Num a) => Vector a -> Vector a -> Vector a
(.+)
{-# INLINE sum #-}
{-# SPECIALIZE INLINE sum :: NonEmpty (Vector Default) -> Vector Default #-}
sort :: (Unbox a, Ord a) => Vector a -> Vector a
sort :: forall a. (Unbox a, Ord a) => Vector a -> Vector a
sort = (forall s. MVector s a -> ST s ()) -> Vector a -> Vector a
forall a.
Unbox a =>
(forall s. MVector s a -> ST s ()) -> Vector a -> Vector a
modify MVector s a -> ST s ()
MVector (PrimState (ST s)) a -> ST s ()
forall s. MVector s a -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) e.
(PrimMonad m, MVector v e, Ord e) =>
v (PrimState m) e -> m ()
Vector.sort
{-# SPECIALIZE sort :: Vector Default -> Vector Default #-}
sortBy :: Unbox a => (a -> a -> Ordering) -> Vector a -> Vector a
sortBy :: forall a. Unbox a => (a -> a -> Ordering) -> Vector a -> Vector a
sortBy a -> a -> Ordering
cmp = (forall s. MVector s a -> ST s ()) -> Vector a -> Vector a
forall a.
Unbox a =>
(forall s. MVector s a -> ST s ()) -> Vector a -> Vector a
modify ((a -> a -> Ordering) -> MVector (PrimState (ST s)) a -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) e.
(PrimMonad m, MVector v e) =>
Comparison e -> v (PrimState m) e -> m ()
Vector.sortBy a -> a -> Ordering
cmp)
sortOn :: (Unbox a, Ord b) => (a -> b) -> Vector a -> Vector a
sortOn :: forall a b. (Unbox a, Ord b) => (a -> b) -> Vector a -> Vector a
sortOn a -> b
key = (a -> a -> Ordering) -> Vector a -> Vector a
forall a. Unbox a => (a -> a -> Ordering) -> Vector a -> Vector a
sortBy (b -> b -> Ordering
forall a. Ord a => a -> a -> Ordering
compare (b -> b -> Ordering) -> (a -> b) -> a -> a -> Ordering
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
`on` a -> b
key)
{-# SPECIALIZE sortOn :: Unbox a => (a -> Default) -> Vector a -> Vector a #-}
argSort :: (Unbox a, Ord a) => Vector a -> Vector Int
argSort :: forall a. (Unbox a, Ord a) => Vector a -> Vector Int
argSort Vector a
vec = (forall s. ST s (MVector s Int)) -> Vector Int
forall a. Unbox a => (forall s. ST s (MVector s a)) -> Vector a
create ((forall s. ST s (MVector s Int)) -> Vector Int)
-> (forall s. ST s (MVector s Int)) -> Vector Int
forall a b. (a -> b) -> a -> b
$ do
MVector s Int
index <- Int -> (Int -> Int) -> ST s (MVector (PrimState (ST s)) Int)
forall (m :: * -> *) a.
(PrimMonad m, Unbox a) =>
Int -> (Int -> a) -> m (MVector (PrimState m) a)
generate (Vector a -> Int
forall a. Unbox a => Vector a -> Int
length Vector a
vec) Int -> Int
forall a. a -> a
id
Comparison Int -> MVector (PrimState (ST s)) Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) e.
(PrimMonad m, MVector v e) =>
Comparison e -> v (PrimState m) e -> m ()
Vector.sortBy (a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare (a -> a -> Ordering) -> (Int -> a) -> Comparison Int
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
`on` Vector a -> Int -> a
forall a. Unbox a => Vector a -> Int -> a
unsafeIndex Vector a
vec) MVector s Int
MVector (PrimState (ST s)) Int
index
pure MVector s Int
index
{-# SPECIALIZE argSort :: Vector Default -> Vector Int #-}
{-# SPECIALIZE argSort :: Vector Int -> Vector Int #-}
sortLike :: (Unbox a, Unbox b, Ord b) => Vector a -> Vector b -> Vector a
sortLike :: forall a b.
(Unbox a, Unbox b, Ord b) =>
Vector a -> Vector b -> Vector a
sortLike = (Vector a -> Vector b -> Vector a)
-> Vector a -> Vector b -> Vector a
forall a b c.
(Unbox a, Unbox b) =>
(Vector a -> Vector b -> c) -> Vector a -> Vector b -> c
withSameLength ((Vector a -> Vector b -> Vector a)
-> Vector a -> Vector b -> Vector a)
-> (Vector a -> Vector b -> Vector a)
-> Vector a
-> Vector b
-> Vector a
forall a b. (a -> b) -> a -> b
$ \Vector a
x Vector b
y ->
Vector a -> Vector Int -> Vector a
forall a. Unbox a => Vector a -> Vector Int -> Vector a
unsafeBackpermute Vector a
x (Vector b -> Vector Int
forall a. (Unbox a, Ord a) => Vector a -> Vector Int
argSort Vector b
y)
{-# SPECIALIZE sortLike :: Unbox a => Vector a -> Vector Default -> Vector a #-}
{-# SPECIALIZE sortLike :: Vector Default -> Vector Int -> Vector Default #-}
normalized :: (Unbox a, Fractional a, Ord a) => Vector a -> Vector a
normalized :: forall a. (Unbox a, Fractional a, Ord a) => Vector a -> Vector a
normalized Vector a
vector
| Vector a -> Bool
forall a. Unbox a => Vector a -> Bool
null Vector a
vector Bool -> Bool -> Bool
|| a
absMax a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
0 = Vector a
vector
| Bool
otherwise = (a -> a) -> Vector a -> Vector a
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
map (a -> a -> a
forall a. Fractional a => a -> a -> a
/ a
absMax) Vector a
vector
where
absMax :: a
absMax = a -> a
forall a. Num a => a -> a
abs ((a -> a) -> Vector a -> a
forall b (v :: * -> *) a.
(Ord b, Vector v a) =>
(a -> b) -> v a -> a
Vector.maximumOn a -> a
forall a. Num a => a -> a
abs Vector a
vector)
{-# SPECIALIZE normalized :: Vector Default -> Vector Default #-}
mean :: (Unbox a, Fractional a) => Vector a -> a
mean :: forall a. (Unbox a, Fractional a) => Vector a -> a
mean Vector a
vector
| Vector a -> Bool
forall a. Unbox a => Vector a -> Bool
null Vector a
vector = a
0
| Bool
otherwise = Vector a -> a
forall (v :: * -> *) a. (Vector v a, Num a) => v a -> a
Vector.sum Vector a
vector a -> a -> a
forall a. Fractional a => a -> a -> a
/ Int -> a
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Vector a -> Int
forall a. Unbox a => Vector a -> Int
length Vector a
vector)
{-# SPECIALIZE mean :: Vector Default -> Default #-}
variance :: (Unbox a, Fractional a) => Vector a -> a
variance :: forall a. (Unbox a, Fractional a) => Vector a -> a
variance Vector a
vector = Vector a -> a
forall a. (Unbox a, Fractional a) => Vector a -> a
mean (Vector a
dev Vector a -> Vector a -> Vector a
forall a. (Unbox a, Num a) => Vector a -> Vector a -> Vector a
.* Vector a
dev)
where
u :: a
u = Vector a -> a
forall a. (Unbox a, Fractional a) => Vector a -> a
mean Vector a
vector
dev :: Vector a
dev = (a -> a) -> Vector a -> Vector a
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
map (\a
x -> a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
u) Vector a
vector
{-# SPECIALIZE variance :: Vector Default -> Default #-}
stdev :: (Unbox a, Floating a) => Vector a -> a
stdev :: forall a. (Unbox a, Floating a) => Vector a -> a
stdev = a -> a
forall a. Floating a => a -> a
sqrt (a -> a) -> (Vector a -> a) -> Vector a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector a -> a
forall a. (Unbox a, Fractional a) => Vector a -> a
variance
{-# SPECIALIZE stdev :: Vector Default -> Default #-}
standardized :: (Unbox a, Floating a, Eq a) => Vector a -> Vector a
standardized :: forall a. (Unbox a, Floating a, Eq a) => Vector a -> Vector a
standardized Vector a
vector
| Vector a -> Bool
forall a. Unbox a => Vector a -> Bool
null Vector a
vector = Vector a
vector
| a
s a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
0 = (a -> a) -> Vector a -> Vector a
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
map (\a
x -> a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
u) Vector a
vector
| Bool
otherwise = (a -> a) -> Vector a -> Vector a
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
map (\a
x -> (a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
u) a -> a -> a
forall a. Fractional a => a -> a -> a
/ a
s) Vector a
vector
where
u :: a
u = Vector a -> a
forall a. (Unbox a, Fractional a) => Vector a -> a
mean Vector a
vector
s :: a
s = Vector a -> a
forall a. (Unbox a, Floating a) => Vector a -> a
stdev Vector a
vector
{-# SPECIALIZE standardized :: Vector Default -> Vector Default #-}
sumM :: (Unbox a, Num a, Monad m) => NonEmpty (m (Vector a)) -> m (Vector a)
sumM :: forall a (m :: * -> *).
(Unbox a, Num a, Monad m) =>
NonEmpty (m (Vector a)) -> m (Vector a)
sumM NonEmpty (m (Vector a))
values = NonEmpty (Vector a) -> Vector a
forall a. (Unbox a, Num a) => NonEmpty (Vector a) -> Vector a
sum (NonEmpty (Vector a) -> Vector a)
-> m (NonEmpty (Vector a)) -> m (Vector a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> NonEmpty (m (Vector a)) -> m (NonEmpty (Vector a))
forall (t :: * -> *) (m :: * -> *) a.
(Traversable t, Monad m) =>
t (m a) -> m (t a)
forall (m :: * -> *) a. Monad m => NonEmpty (m a) -> m (NonEmpty a)
sequence NonEmpty (m (Vector a))
values
{-# INLINE sumM #-}
{-# SPECIALIZE INLINE sumM :: Monad m => NonEmpty (m (Vector Default)) -> m (Vector Default) #-}
{-# SPECIALIZE INLINE sumM :: NonEmpty (Random (Vector Default)) -> Random (Vector Default) #-}
uniformRN :: (Unbox a, Variate a) => a -> a -> Int -> Random (Vector a)
uniformRN :: forall a.
(Unbox a, Variate a) =>
a -> a -> Int -> Random (Vector a)
uniformRN a
lo a
hi Int
count = Int -> Random a -> Random (Vector a)
forall (m :: * -> *) a.
(Monad m, Unbox a) =>
Int -> m a -> m (Vector a)
replicateM Int
count (a -> a -> Random a
forall a. Variate a => a -> a -> Random a
uniformR a
lo a
hi)
{-# SPECIALIZE uniformRN :: Default -> Default -> Int -> Random (Vector Default) #-}
uniformN :: (Unbox a, Variate a, Num a) => Int -> Random (Vector a)
uniformN :: forall a. (Unbox a, Variate a, Num a) => Int -> Random (Vector a)
uniformN = a -> a -> Int -> Random (Vector a)
forall a.
(Unbox a, Variate a) =>
a -> a -> Int -> Random (Vector a)
uniformRN a
0 a
1
{-# SPECIALIZE uniformN :: Int -> Random (Vector Default) #-}
uniformSN :: (Unbox a, Variate a, Num a) => Int -> Random (Vector a)
uniformSN :: forall a. (Unbox a, Variate a, Num a) => Int -> Random (Vector a)
uniformSN = a -> a -> Int -> Random (Vector a)
forall a.
(Unbox a, Variate a) =>
a -> a -> Int -> Random (Vector a)
uniformRN (-a
1) a
1
{-# SPECIALIZE uniformSN :: Int -> Random (Vector Default) #-}
weighted :: (Unbox a, Variate a, Num a) => a -> Vector a -> Random (Vector a)
weighted :: forall a.
(Unbox a, Variate a, Num a) =>
a -> Vector a -> Random (Vector a)
weighted a
maxWeight Vector a
vec = do
a
k <- a -> a -> Random a
forall a. Variate a => a -> a -> Random a
uniformR a
0 a
maxWeight
pure (a
k a -> Vector a -> Vector a
forall a. (Unbox a, Num a) => a -> Vector a -> Vector a
.*. Vector a
vec)
{-# SPECIALIZE weighted :: Default -> Vector Default -> Random (Vector Default) #-}
weightedN :: (Unbox a, Variate a, Num a) => a -> Vector a -> Random (Vector a)
weightedN :: forall a.
(Unbox a, Variate a, Num a) =>
a -> Vector a -> Random (Vector a)
weightedN a
maxWeight Vector a
vec = do
Vector a
k <- a -> a -> Int -> Random (Vector a)
forall a.
(Unbox a, Variate a) =>
a -> a -> Int -> Random (Vector a)
uniformRN a
0 a
maxWeight (Vector a -> Int
forall a. Unbox a => Vector a -> Int
length Vector a
vec)
pure (Vector a
k Vector a -> Vector a -> Vector a
forall a. (Unbox a, Num a) => Vector a -> Vector a -> Vector a
.* Vector a
vec)
{-# SPECIALIZE weightedN :: Default -> Vector Default -> Random (Vector Default) #-}
choice :: NonEmpty (Double, Random a) -> Random a
choice :: forall a. NonEmpty (Double, Random a) -> Random a
choice NonEmpty (Double, Random a)
options = Random (Random a) -> Random a
forall (m :: * -> *) a. Monad m => m (m a) -> m a
join (NonEmpty (Double, Random a) -> Random (Random a)
forall a. NonEmpty (Double, a) -> Random a
weightedChoice NonEmpty (Double, Random a)
options)
{-# INLINE choice #-}