Covariance.hs

module Covariance where

import Means
import qualified Control.Foldl as L
import Control.Applicative

covariance :: (Fractional a, Real a1) => [a1] -> [a1] -> a
covariance xs ys = mean (zipWith (*) xs ys) - mean xs * mean ys

covariance' :: Fractional a => [a] -> [a] -> a
covariance' xs ys = mean' (zipWith (*) xs ys) - mean' xs * mean' ys

covariance'' :: [Double] -> [Double] -> Double
covariance'' xs ys = mean'' (zipWith (*) xs ys) - mean'' xs * mean'' ys

covariance''' :: [Double] -> [Double] -> Double
covariance''' xs ys =
    let mx = mean'' xs
        my = mean'' ys
    in
    sum (zipWith (\x y -> (x - mx) * (y - my)) xs ys) / fromIntegral (length xs)

covariance'''' :: [Double] -> [Double] -> Double
covariance'''' xs ys = mean''' (zipWith (*) xs ys) - mean''' xs * mean''' ys

covariance''''' :: [Double] -> [Double] -> Double
covariance''''' xs ys = let mx = mean''' xs
                            my = mean''' ys
                        in
                        mean''' (zipWith (\x y -> (x - mx) * (y - my)) xs ys)

covariance'''''' :: [Double] -> [Double] -> Double
covariance'''''' xs ys = mean'''' (zipWith (*) xs ys) - mean'''' xs * mean'''' ys

covariance''''''' :: [Double] -> [Double] -> Double
covariance''''''' xs ys = let mx = mean'''' xs
                              my = mean'''' ys
                        in
                        mean'''' (zipWith (\x y -> (x - mx) * (y - my)) xs ys)

Main-weigh.hs

{-# LANGUAGE BangPatterns #-}
module Main where

import Means
import Covariance
import Weigh

main :: IO ()
main = mainWith (do covi
                    covd
                    covFoldl
                    covFolds)

v1 = [1000000..2000000]

v2 = [2000000..3000000]

covi :: Weigh ()
covi = do func "naive mean" cov1 v2
          func "applicative mean" cov2 v2
       where cov1 = covariance v1 :: [Double] -> Double
             cov2 = covariance' v1 :: [Double] -> Double

covd :: Weigh ()
covd = do func "optimized mean, naive zipWith" cov3 v2
          func "optimized mean, hand-tuned zipWith" cov4 v2
       where cov3 = covariance'' v1 :: [Double] -> Double
             cov4 = covariance''' v1 :: [Double] -> Double

covFoldl :: Weigh ()
covFoldl = do func "foldl mean" cov5 v2
              func "foldl mean, tuned zipWith" cov6 v2
           where cov5 = covariance'''' v1 :: [Double] -> Double
                 cov6 = covariance''''' v1 :: [Double] -> Double

covFolds :: Weigh ()
covFolds = do func "folds mean" cov7 v2
              func "folds mean, tuned zipWith" cov8 v2
           where cov7 = covariance'''''' v1 :: [Double] -> Double
                 cov8 = covariance''''''' v1 :: [Double] -> Double


{-
naive mean                          723,716,168  1,382
applicative mean                    723,714,736  1,382
optimized mean, naive zipWith       456,000,688    875
optimized mean, hand-tuned zipWith  336,000,592    642
foldl mean                          336,000,568    642
foldl mean, tuned zipWith           336,000,568    642
folds mean                          456,000,784    875
folds mean, tuned zipWith           456,000,888    871
-}

Main.hs

{-# LANGUAGE BangPatterns #-}
module Main where

import Means
import Covariance
import Weigh

main :: IO ()
main = do
         --let cov = covariance v1 v2
         --let cov = covariance' v1 v2
         --let cov = covariance'' v1 v2
         --let cov = covariance''' v1 v2
         let cov = covariance'''' v1 v2
         --let cov = covariance''''' v1 v2
         --let cov = covariance'''''' v1 v2
         --let cov = covariance''''''' v1 v2
         print cov

v1 = [1000000..2000000]

v2 = [2000000..3000000]

Means.hs

module Means where

import Data.List
import Control.Applicative
import qualified Control.Foldl as L
import qualified Data.Fold as F

mean :: (Fractional a1, Real a, Foldable t) => t a -> a1
mean xs = realToFrac (sum xs) / realToFrac (length xs)

mean' :: (Foldable f, Fractional a) => f a -> a
mean' = liftA2 (/) sum (fromIntegral . length)

data Pair = Pair {-# UNPACK #-}!Int {-# UNPACK #-}!Double

mean'' :: [Double] -> Double
mean'' xs = s / fromIntegral n
  where
    Pair n s = foldl' k (Pair 0 0) xs
    k (Pair n s) x = Pair (n+1) (s+x)

mean''' :: [Double] -> Double
mean''' = L.fold (liftA2 (/) L.sum L.genericLength)

sumL :: F.L Double Double
sumL = F.L id (+) 0

lengthL :: F.L a Double
lengthL = F.L id (const . (+1)) 0

mean'''' :: [Double] -> Double
mean'''' = flip F.run (liftA2 (/) sumL lengthL)