hstatistics 0.2.5.2 → 0.2.5.3

raw patch · 9 files changed

+111/−84 lines, 9 filesdep ~basedep ~hmatrixdep ~hmatrix-gsl-statsPVP: major bump suggested

API removals or changes: PVP suggests a major version bump

Dependency ranges changed: base, hmatrix, hmatrix-gsl-stats

API changes (from Hackage documentation)

Files

• CHANGES +3/−0
• hstatistics.cabal +6/−5
• lib/Numeric/Statistics.hs +32/−32
• lib/Numeric/Statistics/Histogram.hs +9/−8
• lib/Numeric/Statistics/ICA.hs +27/−12
• lib/Numeric/Statistics/Information.hs +14/−8
• lib/Numeric/Statistics/PCA.hs +8/−8
• lib/Numeric/Statistics/PDF.hs +4/−3
• lib/Numeric/Statistics/Surrogate.hs +8/−8

@@ -103,3 +103,6 @@  0.2.5.2: 		exposed pcaReduceN as requested by Tom Nielsen++0.2.5.3:+		updated for hmatrix 0.17

@@ -1,8 +1,8 @@ Name:               hstatistics-Version:            0.2.5.2+Version:            0.2.5.3 License:            BSD3 License-file:       LICENSE-Copyright:          (c) A.V.H. McPhail 2010, 2011, 2012, 2013+Copyright:          (c) A.V.H. McPhail 2010, 2011, 2012, 2013, 2014 Author:             Vivian McPhail Maintainer:         haskell.vivian.mcphail <at> gmail <dot> com Stability:          provisional@@ -16,7 +16,7 @@      .      Feature requests, suggestions, and bug fixes welcome. Category:           Math, Statistics-tested-with:        GHC ==7.6.3+tested-with:        GHC ==7.10.2  cabal-version:      >=1.8 @@ -30,8 +30,8 @@     Build-Depends:      base >= 4 && < 5,                         array, random,                         vector,-                        hmatrix >= 0.10.0.0,-                        hmatrix-gsl-stats >= 0.1.2.9+                        hmatrix >= 0.17,+                        hmatrix-gsl-stats >= 0.4      Extensions:          @@ -51,6 +51,7 @@     ghc-options:        -Wall -fno-warn-missing-signatures                               -fno-warn-orphans                               -fno-warn-unused-binds+                        -O2  source-repository head     type:     git

@@ -2,8 +2,8 @@ ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics--- Copyright   :  (c) A. V. H. McPhail 2010, 2012--- License     :  BSD+-- Copyright   :  (c) A. V. H. McPhail 2010, 2012, 2014+-- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com -- Stability   :  provisional@@ -31,14 +31,14 @@  ----------------------------------------------------------------------------- ---import Numeric.Vector---import Numeric.Matrix---import Numeric.Container-import Numeric.LinearAlgebra+import Numeric.LinearAlgebra hiding(range)+--import Numeric.LinearAlgebra.Data hiding(range)+--import Numeric.LinearAlgebra.Devel  import qualified Data.Array.IArray as I  import qualified Data.List as DL import qualified Data.Vector.Generic as GV+--import qualified Data.Vector.Storable as SV  import Foreign.Storable @@ -68,36 +68,36 @@ -----------------------------------------------------------------------------  -- | the mean of a list of vectors-meanList :: (Container Vector a, Num (Vector a)) => [Sample a] -> Sample a+meanList :: (Container Vector a, Num (Vector a), Fractional a) => [Sample a] -> Sample a meanList []     = error "meanVectors: empty list" meanList [s]    = s meanList (s:ss) = let ln = fromIntegral $ length ss + 1                   in scale (recip ln) $ foldl (+) s ss  -- | the mean of an array of vectors-meanArray :: (Container Vector a, Num (Vector a)) => Samples a -> Sample a+meanArray :: (Container Vector a, Num (Vector a), Fractional a) => Samples a -> Sample a meanArray a = meanList $ I.elems a  -- | the mean of a matrix with data series in rows-meanMatrix :: (Container Vector a, Num (Vector a), Element a) => Matrix a -> Sample a+meanMatrix :: (Container Vector a, Num (Vector a), Element a, Fractional a) => Matrix a -> Sample a meanMatrix a = meanList $ toRows a  -----------------------------------------------------------------------------  -- | the variance of a list of vectors-varianceList :: (Container Vector a, Floating (Vector a)) => [Sample a] -> Sample a+varianceList :: (Container Vector a, Floating (Vector a), Num a, Fractional a) => [Sample a] -> Sample a varianceList []  = error "varianceList: empty list"-varianceList [s] = constant 0 (dim s)+varianceList [s] = konst 0 (size s) varianceList l   = let mxs = meanList (map (** 2) l)                        msx = (meanList l) ** 2                    in mxs - msx  -- | the variance of an array of vectors-varianceArray :: (Container Vector a, Floating (Vector a)) => Samples a -> Sample a+varianceArray :: (Container Vector a, Floating (Vector a), Num a, Fractional a) => Samples a -> Sample a varianceArray a = varianceList $ I.elems a  -- | the variance of a matrix with data series in rows-varianceMatrix :: (Container Vector a, Floating (Vector a), Element a) => Matrix a -> Sample a+varianceMatrix :: (Container Vector a, Floating (Vector a), Element a, Num a, Fractional a) => Matrix a -> Sample a varianceMatrix a = varianceList $ toRows a  -----------------------------------------------------------------------------@@ -128,11 +128,11 @@     -> Vector Double -- ^ intervals     -> Vector Int    -- ^ data indexed by bin cut v c  = let c' = sort c-           in mapVector (\x -> cut_helper 0 x c') v +           in GV.map (\x -> cut_helper 0 x c') v      where       cut_helper j x d -          | j >= dim d                       = error "Numeric.Statistics: cut: data point not within interval"-          | x >= (d @> j) && x <= (d @> (j+1)) = j+          | j >= size d                       = error "Numeric.Statistics: cut: data point not within interval"+          | x >= (d `atIndex` j) && x <= (d `atIndex` (j+1)) = j           | otherwise                       = cut_helper (j + 1) x d  -----------------------------------------------------------------------------@@ -142,19 +142,19 @@ --ranks :: Vector Double -> Vector Double ranks :: (Fractional b, Storable b) => Vector Double -> Vector b ranks v = let v' = sort v-          in mapVector (\x -> 1 + rank_helper x v') v+          in GV.map (\x -> 1 + rank_helper x v') v               where rank_helper x v' = let is = GV.elemIndices x v'-                                       in (realToFrac (GV.foldl (+) 0 is)) / (fromIntegral $ dim is)+                                       in (realToFrac (GV.foldl (+) 0 is)) / (fromIntegral $ GV.length is)  -----------------------------------------------------------------------------  -- | kendall's rank correlation τ kendall :: Vector Double -> Vector Double -> Matrix Double-kendall x y = let ln = dim x+kendall x y = let ln = size x                   rx = ranks x                   ry = ranks y                   r = fromColumns [rx,ry]-                  m = signum $ (kronecker r (asColumn $ constant 1.0 ln)) - (kronecker (asRow $ constant 1.0 ln) r)+                  m = signum $ (kronecker r (asColumn $ konst 1.0 ln)) - (kronecker (asRow $ konst 1.0 ln) r)                   c = rows m - 1               in correlationCoefficientMatrix $ I.listArray (0,c) (toColumns m) @@ -164,7 +164,7 @@ --logit :: Vector Double -> Vector Double logit :: (Floating b, Storable b)         => Vector b -> Vector b-logit v =  mapVector (\x -> - (log ((1 / x) - 1))) v+logit v =  GV.map (\x -> - (log ((1 / x) - 1))) v  ----------------------------------------------------------------------------- @@ -181,16 +181,16 @@                                  Just m  -> m                       xm     = fromRows $ map ((-) xu) $ toRows $ fromColumns xl                       --um     = asColumn xu-                      --w      = ((trans xm) <> xm + (trans um) <> um)/(fromIntegral $ xr - 1)+                      --w      = ((tr' xm) <> xm + (tr' um) <> um)/(fromIntegral $ xr - 1)                       --w'     = inv w-                  in ((xm <> s' <> (trans xm)) @@> (0,0)) +                  in ((xm <> s' <> (tr' xm)) `atIndex` (0,0))   -----------------------------------------------------------------------------  -- | a list of element frequencies mode :: Vector Double -> [(Double,Integer)] mode v = let w = sort v-         in DL.sortBy (\(_,n) (_,n') -> compare n' n) $ foldVector freqs [] w+         in DL.sortBy (\(_,n) (_,n') -> compare n' n) $ GV.foldr freqs [] w             where freqs x []          = [(x,1)]                   freqs x ((f,n):fns)                       | f == x         = ((f,n+1):fns) @@ -211,7 +211,7 @@ --    | p == 2     = variance v -- gives sample variance     | otherwise = let u = if c then centre v else v                       w = if a then abs u else u-                      x = mapVector (** (fromIntegral p)) w+                      x = GV.map (** (fromIntegral p)) w                   in mean x  -----------------------------------------------------------------------------@@ -227,13 +227,13 @@     | rows x /= rows y = error "Numeric.Statistics: ols: incorrect matrix dimensions"     | otherwise       = let (xr,xc) = (rows x,cols x)                             (yr,yc) = (rows y,cols y)-                            z = (trans x) <> x+                            z = (tr' x) <> x                             r = rank z                             beta = if r == xc -                                      then (inv z) <> (trans x) <> y+                                      then (inv z) <> (tr' x) <> y                                       else (pinv x) <> y                             rr = y - x <> beta-                            sigma = ((trans rr) <> rr) / (fromIntegral $ xr - r)+                            sigma = ((tr' rr) <> rr) / (fromIntegral $ xr - r)                         in (beta,rr,sigma)  -----------------------------------------------------------------------------@@ -247,18 +247,18 @@ -----------------------------------------------------------------------------  -- | the difference between the maximum and minimum of the input-range :: Container c e => c e -> e+range :: (Container c e, Num e) => c e -> e range v = maxElement v - minElement v  -----------------------------------------------------------------------------  -- | count the number of runs greater than or equal to @n@ in the data-run_count :: (Num a, Num t, Ord b, Ord a, Storable b) +run_count :: (Num a, Num t, Ord b, Ord a, Storable b, Container Vector b)              => a             -- ^ longest run to count           -> Vector b        -- ^ data           -> [(a, t)]        -- ^ [(run length,count)]-run_count n v = let w = subVector 1 (dim v - 1) v-                    x = foldVector run_count' [(1,v @> 0)] w+run_count n v = let w = subVector 1 (size v - 1) v+                    x = GV.foldr run_count' [(1,v `atIndex` 0)] w                     y = map fst x                     z = takeWhile (<= n) $  DL.sort y                 in foldr count [] z

@@ -1,8 +1,8 @@ ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics.Histogram--- Copyright   :  (c) A. V. H. McPhail 2010--- License     :  GPL-style+-- Copyright   :  (c) A. V. H. McPhail 2010, 2014+-- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com -- Stability   :  provisional@@ -20,7 +20,8 @@  ----------------------------------------------------------------------------- -import Data.Packed.Vector+--import qualified Data.Packed.Vector as V+import qualified Data.Vector.Storable as V  import qualified Numeric.GSL.Histogram as H --import qualified Numeric.GSL.Histogram2D as H2@@ -32,7 +33,7 @@  ----------------------------------------------------------------------------- -vectorToTuples = toTuples . toList+vectorToTuples = toTuples . V.toList     where toTuples []         = error "need a minimum of two elements"           toTuples [_]        = error "need a minimum of two elements"           toTuples [x1,x2]    = [(x1,x2)]@@ -40,14 +41,14 @@  ----------------------------------------------------------------------------- -cumulativeToHistogram :: (Double -> Double)  -- ^ the cumulative distribution function D(x <= X)-                   -> Vector Double          -- ^ the bins+cumulativeToHistogram :: (Double -> Double)   -- ^ the cumulative distribution function D(x <= X)+                   -> V.Vector Double        -- ^ the bins                    -> H.Histogram            -- ^ the resulting histogram cumulativeToHistogram f v = H.addListWeighted (H.emptyRanges v) $ map (\(x1,x2) -> ((x1 + x2) / 2.0,f x2 - f x1)) (vectorToTuples v) -gaussianHistogram :: Double              -- ^ mean+gaussianHistogram :: Double               -- ^ mean                   -> Double              -- ^ standard deviation-                  -> Vector Double       -- ^ the bins+                  -> V.Vector Double     -- ^ the bins                   -> H.Histogram         -- ^ the resulting histogram gaussianHistogram u s = cumulativeToHistogram (\x -> C.density_1p C.Gaussian C.Lower s (x-u)) 

@@ -1,7 +1,8 @@+{-# LANGUAGE FlexibleContexts #-} ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics.ICA--- Copyright   :  (c) A. V. H. McPhail 2010+-- Copyright   :  (c) A. V. H. McPhail 2010, 2014 -- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com@@ -23,7 +24,9 @@ module Numeric.Statistics.ICA (                                sigmoid, sigmoid',                                demean, whiten,-                               ica, icaDefaults+                               ica, icaDefaults,+                               --+                               NormType(..)                           ) where  @@ -33,6 +36,8 @@  import Numeric.LinearAlgebra +import qualified Data.Vector.Generic as GV+ import Numeric.GSL.Statistics  import Numeric.Statistics@@ -41,6 +46,16 @@  ----------------------------------------------------------------------------- +data NormType = NormZero | NormOne | NormTwo | NormInf++pnorm :: Normed (Vector a) => NormType -> Vector a -> R+pnorm NormZero = norm_0+pnorm NormOne  = norm_1+pnorm NormTwo  = norm_2+pnorm NormInf  = norm_Inf++-----------------------------------------------------------------------------+ -- | sigmoid transfer function sigmoid :: Double -> Double sigmoid u = u * exp((-u**2)/2)@@ -76,7 +91,7 @@        -> Double                                      -- ^ eigenvalue threshold        -> (I.Array Int (Vector Double),Matrix Double) -- ^ (whitened data,transform) whiten d q = let cv = covarianceMatrix d-                 (val',vec') = eigSH cv           -- the covariance matrix is real symmetric+                 (val',vec') = eigSH $ trustSym cv  -- the covariance matrix is real symmetric                  val = toList val'                  vec = toColumns vec'                  v' = zip val vec@@ -85,7 +100,7 @@                  dd = diag $ (** (-0.5)) $ fromList dd'  -- square root of eigenvalues diagonalised                  e = fromColumns e'                  x = fromRows $ I.elems d-                 t = e <> dd <> trans e          -- the actual mathematics+                 t = e <> dd <> tr' e          -- the actual mathematics                  x' = t <> x                     -- the actual mathematics                  d' = I.listArray (I.bounds d) (toRows x')               in (d',t)@@ -125,23 +140,23 @@ update :: (Double -> Double) -> (Double -> Double) -> Matrix Double -> Matrix Double -> Matrix Double update g g' w x = let y = w <> x                       ys = toRows y-                      bis = map (\y' -> - mean (y' * (mapVector g y'))) ys-                      ais = zipWith (\b y' -> -1 / (b - mean (mapVector g y'))) bis ys+                      bis = map (\y' -> - mean (y' * (GV.map g y'))) ys+                      ais = zipWith (\b y' -> -1 / (b - mean (GV.map g y'))) bis ys                       r = rows y                       ix = ((1,1),(r,r))-                      cov = fromArray2D $ I.listArray ix $ map (\(m,n) -> covariance (mapVector g' (ys!!(m-1))) (ys!!(n-1))) $ I.range ix+                      cov = fromArray2D $ I.listArray ix $ map (\(m,n) -> covariance (GV.map g' (ys!!(m-1))) (ys!!(n-1))) $ I.range ix                   in w + (diag $ fromList ais) <> ((diag $ fromList bis) + cov) <> w    decorrelate :: Matrix Double -> Matrix Double decorrelate m = let (d',v') = eig m                     d = fst $ fromComplex d'                     v = fst $ fromComplex v'-                in v <> (diag (d ** (-0.5))) <> trans v <> m-{-decorrelate n t w = let w' = w / (scalar $ sqrt $ pnorm n (w <> trans w))+                in v <> (diag (d ** (-0.5))) <> tr' v <> m+{-decorrelate n t w = let w' = w / (scalar $ sqrt $ pnorm n (w <> tr' w))                     in decorrelate' t w w'     where decorrelate' t' m m'                | converged t' m m' = m'-              | otherwise         = decorrelate' t' m' ((scale 1.5 m') - (scale 0.5 (m' <> trans m' <> m')))+              | otherwise         = decorrelate' t' m' ((scale 1.5 m') - (scale 0.5 (m' <> tr' m' <> m'))) -}  normalise :: NormType -> Matrix Double -> Matrix Double@@ -192,7 +207,7 @@             -> I.Array Int (Vector Double) -- ^ data             -> (I.Array Int (Vector Double),Matrix Double) -- ^ transformed data, ica transform icaDefaults r a = let c = I.rangeSize $ I.bounds a-                      s = (dim $ (a I.! 1)) `div` 16-                  in ica r sigmoid sigmoid' Infinity 0.0000001 s a+                      s = (GV.length $ (a I.! 1)) `div` 16+                  in ica r sigmoid sigmoid' NormInf 0.0000001 s a  -----------------------------------------------------------------------------

@@ -2,7 +2,7 @@ ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics.Information--- Copyright   :  (c) A. V. H. McPhail 2010+-- Copyright   :  (c) A. V. H. McPhail 2010, 2014 -- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com@@ -23,15 +23,21 @@  import Numeric.Statistics.PDF -import Numeric.LinearAlgebra+import qualified Numeric.LinearAlgebra as LA+--import Numeric.LinearAlgebra.Data hiding(Vector) +--import qualified Data.Vector as DV+import qualified Data.Vector.Storable as V++import Prelude hiding(map,zip)+ -----------------------------------------------------------------------------  zeroToOne x     | x == 0.0  = 1.0     | otherwise = x -logE = mapVector (log . zeroToOne)+logE = V.map (log . zeroToOne)   -----------------------------------------------------------------------------@@ -39,21 +45,21 @@ -- | the entropy \sum p_i l\ln{p_i} of a sequence entropy :: PDF a Double          => a                       -- ^ the underlying distribution-        -> Vector Double           -- ^ the sequence+        -> LA.Vector Double           -- ^ the sequence         -> Double                  -- ^ the entropy entropy p x = let ps = probability p x-              in negate $ (dot ps (logE ps))+              in negate $ (LA.dot ps (logE ps))  -- | the mutual information \sum_x \sum_y p(x,y) \ln{\frac{p(x,y)}{p(x)p(y)}} mutual_information :: (PDF a Double, PDF b (Double,Double))                     => b                                          -- ^ the underlying distribution                    -> a                                          -- ^ the first dimension distribution                    -> a                                          -- ^ the second dimension distribution-                   -> (Vector Double, Vector Double)             -- ^ the sequence+                   -> (LA.Vector Double, LA.Vector Double)             -- ^ the sequence                    -> Double         -- ^ the mutual information-mutual_information p px py (x,y) = let ps = probability p $ zipVector x y+mutual_information p px py (x,y) = let ps = probability p $ V.zipWith (,) x y                                        xs = probability px x                                        ys = probability py y-                                   in (dot ps (logE ps - logE (xs*ys)))+                                   in (LA.dot ps (logE ps - logE (xs*ys)))  -----------------------------------------------------------------------------

@@ -1,7 +1,7 @@ ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics.PCA--- Copyright   :  (c) A. V. H. McPhail 2010+-- Copyright   :  (c) A. V. H. McPhail 2010, 2014 -- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com@@ -38,7 +38,7 @@     -> Matrix Double pca d q = let d' = fmap (\x -> x - (scalar $ mean x)) d -- remove the mean from each dimension               cv = covarianceMatrix d'-              (val',vec') = eigSH cv           -- the covariance matrix is real symmetric+              (val',vec') = eigSH $ trustSym cv -- the covariance matrix is real symmetric               val = toList val'               vec = toColumns vec'               v' = zip val vec@@ -52,7 +52,7 @@      -> Matrix Double pcaN d n = let d' = fmap (\x -> x - (scalar $ mean x)) d -- remove the mean from each dimension                cv = covarianceMatrix d'-               (val',vec') = eigSH cv           -- the covariance matrix is real symmetric+               (val',vec') = eigSH $ trustSym cv  -- the covariance matrix is real symmetric                val = toList val'                vec = toColumns vec'                v' = zip val vec@@ -65,7 +65,7 @@              -> Matrix Double                  -- ^ the principal components              -> I.Array Int (Vector Double)    -- ^ the transformed data pcaTransform d m = let d' = fmap (\x -> x - (scalar $ mean x)) d -- remove the mean from each dimension-                   in I.listArray (1,cols m) $ toRows $ (trans m) <> (fromRows $ I.elems d')+                   in I.listArray (1,cols m) $ toRows $ (tr' m) <> (fromRows $ I.elems d')  -- | perform a dimension-reducing PCA modification,  --     using an eigenvalue threshhold@@ -75,13 +75,13 @@ pcaReduce d q = let u = fmap (scalar . mean) d                     d' = zipWith (-) (I.elems d) (I.elems u)                     cv = covarianceMatrix $ I.listArray (I.bounds d) d'-                    (val',vec') = eigSH cv           -- the covariance matrix is real symmetric+                    (val',vec') = eigSH $ trustSym cv -- the covariance matrix is real symmetric                     val = toList val'                     vec = toColumns vec'                     v' = zip val vec                     v = filter (\(x,_) -> x > q) v'  -- keep only eigens > than parameter                     m = fromColumns $ snd $ unzip v-                 in I.listArray (I.bounds d) $ zipWith (+) (toRows $ m <> (trans m) <> fromRows d') (I.elems u) +                 in I.listArray (I.bounds d) $ zipWith (+) (toRows $ m <> (tr' m) <> fromRows d') (I.elems u)   -- | perform a dimension-reducing PCA modification, using N components pcaReduceN :: I.Array Int (Vector Double)      -- ^ the data@@ -90,12 +90,12 @@ pcaReduceN d n = let u = fmap (scalar . mean) d                      d' = zipWith (-) (I.elems d) (I.elems u)                      cv = covarianceMatrix $ I.listArray (I.bounds d) d'-                     (val',vec') = eigSH cv           -- the covariance matrix is real symmetric+                     (val',vec') = eigSH $ trustSym cv -- the covariance matrix is real symmetric                      val = toList val'                      vec = toColumns vec'                      v' = zip val vec                      v = take n $ reverse $ sortBy (comparing fst) v'                      m = fromColumns $ snd $ unzip v-                  in I.listArray (I.bounds d) $ zipWith (+) (toRows $ m <> (trans m) <> fromRows d') (I.elems u) +                  in I.listArray (I.bounds d) $ zipWith (+) (toRows $ m <> (tr' m) <> fromRows d') (I.elems u)   -----------------------------------------------------------------------------

@@ -3,7 +3,7 @@ ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics.PDF--- Copyright   :  (c) A. V. H. McPhail 2010+-- Copyright   :  (c) A. V. H. McPhail 2010, 2014 -- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com@@ -25,7 +25,8 @@  ----------------------------------------------------------------------------- -import qualified Data.Packed.Vector as V+--import qualified Data.Packed.Vector as V+import qualified Data.Vector.Storable as V  import qualified Numeric.GSL.Histogram as H import qualified Numeric.GSL.Histogram2D as H2@@ -46,7 +47,7 @@     probability :: b -> V.Vector a -> V.Vector Double  instance Storable b => PDF (PDFFunction b) b where-    probability (P_Func f) = V.mapVector f+    probability (P_Func f) = V.map f  instance PDF H.Histogram Double where     probability = H.prob

@@ -1,7 +1,7 @@ ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.Statistics.Surrogate--- Copyright   :  (c) Alexander Vivian Hugh McPhail 2010+-- Copyright   :  (c) Alexander Vivian Hugh McPhail 2010, 2014 -- License     :  BSD3 -- -- Maintainer  :  haskell.vivian.mcphail <at> gmail <dot> com@@ -19,8 +19,8 @@  ----------------------------------------------------------------------------- -import Data.Packed.Vector---import Data.Packed.Matrix+--import qualified Data.Packed.Vector as V+import qualified Data.Vector.Storable as V  import qualified Data.Array.IArray as I  @@ -33,13 +33,13 @@ -- | perform an analysis using surrogate data surrogate :: Int                          -- ^ random seed           -> Int                          -- ^ number of repetitions-          -> (I.Array Int (Vector Double) -> a) -- ^ the evaluation function-          -> I.Array Int (Vector Double)  -- ^ the data+          -> (I.Array Int (V.Vector Double) -> a) -- ^ the evaluation function+          -> I.Array Int (V.Vector Double)  -- ^ the data           -> I.Array Int a                -- ^ the results, with the evaluated real data in position 1 and the rest of the array containing the evaluated surrogate data                                                                                                                surrogate r n f d = I.listArray (1,n+1) $ (f d) : (surrogate' (mkStdGen r) n f d) -surrogate' :: StdGen -> Int-> (I.Array Int (Vector Double) -> a) -> I.Array Int (Vector Double) -> [a]+surrogate' :: StdGen -> Int-> (I.Array Int (V.Vector Double) -> a) -> I.Array Int (V.Vector Double) -> [a] surrogate' _ 0 _ _ = [] surrogate' g n f d = let (g',g'') = split g                          d' = permute_data g' d@@ -50,10 +50,10 @@ randomList g n = let (r,g') = random g                  in r : (randomList g' (n-1)) -permute_data :: StdGen -> I.Array Int (Vector Double) -> I.Array Int (Vector Double)+permute_data :: StdGen -> I.Array Int (V.Vector Double) -> I.Array Int (V.Vector Double) permute_data g d = let s = I.rangeSize $ I.bounds d                        rs = randomList g s                        ds = zip rs $ I.elems d-                   in I.listArray (I.bounds d) $ map (\(r,v) -> permute (random_permute r (dim v)) v) ds+                   in I.listArray (I.bounds d) $ map (\(r,v) -> permute (random_permute r (V.length v)) v) ds   -----------------------------------------------------------------------------