diff --git a/CHANGELOG.md b/CHANGELOG.md new file mode 100644 --- /dev/null +++ b/CHANGELOG.md @@ -0,0 +1,3 @@ +#### 0.0.0 + +* Initial release. diff --git a/LICENSE b/LICENSE new file mode 100644 --- /dev/null +++ b/LICENSE @@ -0,0 +1,26 @@ +Copyright (C) 2018-2019 Stevan Andjelkovic, Robert Danitz + +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + + * Redistributions in binary form must reproduce the above + copyright notice, this list of conditions and the following + disclaimer in the documentation and/or other materials provided + with the distribution. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT +OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, +DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY +THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT +(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. diff --git a/README.md b/README.md new file mode 100644 --- /dev/null +++ b/README.md @@ -0,0 +1,72 @@ +## markov-chain-usage-model + +[![Hackage](https://img.shields.io/hackage/v/markov-chain-usage-model.svg)](https://hackage.haskell.org/package/markov-chain-usage-model) +[![Stackage Nightly](http://stackage.org/package/markov-chain-usage-model/badge/nightly)](http://stackage.org/nightly/package/markov-chain-usage-model) +[![Stackage LTS](http://stackage.org/package/markov-chain-usage-model/badge/lts)](http://stackage.org/lts/package/markov-chain-usage-model) +[![Build Status](https://api.travis-ci.org/advancedtelematic/markov-chain-usage-model.svg?branch=master)](https://travis-ci.org/advancedtelematic/markov-chain-usage-model) + +`markov-chain-usage-model` is a Haskell library for performing computations on +Markov chain based usage models. + +### References + + * [A theory of software reliability and its + application](http://doi.ieeecomputersociety.org/10.1109/TSE.1975.6312856) + (1975) by J. D. Musa; + + * [Engineering software under statistical quality + control](https://doi.org/10.1109/52.60601) (1990) by R. H. Cobb and H. D. + Mills; + + * [Operational profiles in software-reliability + engineering](https://doi.org/10.1109/52.199724) (1993) by J. D. Musa; + + * [A Markov chain model for statistical software + testing](https://doi.org/10.1109/32.328991) (1994) by J. A. Whittaker and M. + G. Thomason; + + * [Statistical testing of software based on a usage + model](https://doi.org/10.1002/spe.4380250106) (1995) by G. H. Walton, J. H. + Poore and C. J. Trammell; + + * [Quantifying the reliability of software: statistical testing based on a + usage model](https://doi.org/10.1109/SESS.1995.525966) (1995) by C. + Trammell; + + * [Computations for Markov Chain Usage + Models](http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.416.2257&rep=rep1&type=pdf) + (2000) by S. J. Prowell; + + * [Computing system reliability using Markov chain usage + models](https://doi.org/10.1016/S0164-1212(03)00241-3) (2004) by S. J. + Prowell and J. H. Poore; + + * Software Reliability Engineering: More Reliable Software Faster and Cheaper + (2004) by J. D. Musa; + + * [A Simpler and More Direct Derivation of System Reliability Using Markov + Chain Usage + Models](https://ksiresearchorg.ipage.com/seke/seke17paper/seke17paper_91.pdf) + (2017) by L. Lin, Y. Xue and F. Song; + + * [On A Simpler and Faster Derivation of Single Use Reliability Mean and + Variance for Model-Based Statistical + Testing](http://ksiresearchorg.ipage.com/seke/seke18paper/seke18paper_26.pdf) + (2018) by Y. Xue, L. Lin, X. Sun and F. Song; + + * [Stopping criteria for statistical + testing](https://doi.org/10.1016/S0950-5849(00)00110-5) (2000) by K. Sayre + and J. H. Poore; + + * [A cost-benefit stopping criterion for statistical + testing](https://doi.org/10.1109/HICSS.2004.1265715) (2004) by S. J. + Prowell. + +### See also + + * The JUMBL library and tool ([documentation](http://jumbl.sourceforge.net/) + and [code](https://sourceforge.net/p/jumbl/code/ci/master/tree/)) + +### License + +BSD-style (see the file LICENSE). diff --git a/Setup.hs b/Setup.hs new file mode 100644 --- /dev/null +++ b/Setup.hs @@ -0,0 +1,2 @@ +import Distribution.Simple +main = defaultMain diff --git a/markov-chain-usage-model.cabal b/markov-chain-usage-model.cabal new file mode 100644 --- /dev/null +++ b/markov-chain-usage-model.cabal @@ -0,0 +1,76 @@ +cabal-version: 1.12 + +name: markov-chain-usage-model +version: 0.0.0 +synopsis: Computations for Markov chain usage models +description: Please see the README on GitHub at +homepage: https://github.com/advancedtelematic/markov-chain-usage-model#readme +bug-reports: https://github.com/advancedtelematic/markov-chain-usage-model/issues +author: Stevan Andjelkovic +maintainer: stevan.andjelkovic@here.com +copyright: Copyright (C) 2018-2019, HERE Europe B.V. +category: Testing +license: BSD2 +license-file: LICENSE +build-type: Simple +tested-with: GHC == 8.2.2, GHC == 8.4.3, GHC == 8.6.3 +extra-source-files: + README.md + CHANGELOG.md + +source-repository head + type: git + location: https://github.com/advancedtelematic/markov-chain-usage-model + +library + exposed-modules: + MarkovChain + other-modules: + Paths_markov_chain_usage_model + hs-source-dirs: + src + ghc-options: + -Weverything + -Wno-missing-exported-signatures + -Wno-missing-import-lists + -Wno-missed-specialisations + -Wno-all-missed-specialisations + -Wno-unsafe + -Wno-safe + -Wno-missing-local-signatures + -Wno-monomorphism-restriction + build-depends: + base >=4.10 && <5 + , matrix >= 0.3.6.1 + , vector >= 0.10 + default-language: Haskell2010 + +test-suite markov-chain-usage-model-test + type: exitcode-stdio-1.0 + main-is: Spec.hs + other-modules: + Paths_markov_chain_usage_model + , Unit + hs-source-dirs: + test + ghc-options: + -threaded -rtsopts -with-rtsopts=-N + -Weverything + -Wno-missing-exported-signatures + -Wno-missing-import-lists + -Wno-missed-specialisations + -Wno-all-missed-specialisations + -Wno-unsafe + -Wno-safe + -Wno-missing-local-signatures + -Wno-monomorphism-restriction + build-depends: + base >=4.7 && <5 + , doctest + , markov-chain-usage-model + , matrix + , tasty + , tasty-discover + , tasty-hunit + , vector + default-language: Haskell2010 diff --git a/src/MarkovChain.hs b/src/MarkovChain.hs new file mode 100644 --- /dev/null +++ b/src/MarkovChain.hs @@ -0,0 +1,418 @@ +module MarkovChain + ( M + , V + , (.*) + , (./) + , (.-) + , (^*) + , norm_1 + , norm_F + , fundamental + , occVariance + , perron + , sigma + , nodeProbabilities + , expectedLength + , expectedArcReliability + , transientReliability + , singleUseReliability + -- , singleUseReliabilityIO + -- , loadStaticMatrix + -- , load2StaticMatrices + -- , saveStaticMatrix + , kullbackLeibler + ) + where + +import Data.Bifunctor + (bimap) +import qualified Data.List as List +import Data.Matrix as Matrix +import Data.Vector + (Vector) +import qualified Data.Vector as Vector +import Prelude hiding + (pi) + +------------------------------------------------------------------------ + +type M = Matrix Double +type V = Vector Double + +-- point-wise +(.*) :: M -> M -> M +(.*) = elementwise (*) +infixl 7 .* + +(./) :: M -> M -> M +(./) = elementwise (/) +infixl 7 ./ + +(.-) :: M -> M -> M +(.-) = elementwise (-) +infixl 6 .- + +-- scalar +(^*) :: Double -> M -> M +(^*) x = fmap (x *) +infixl 7 ^* + +(^/) :: Double -> M -> M +(^/) x = fmap (x /) +infixl 7 ^/ + +norm_1 :: M -> Double +norm_1 = sum . toList + +norm_F :: M -> Double +norm_F = sqrt . sum . toList . fmap (^(2 :: Int)) + +{- $setup + >>> :set -XFlexibleContexts + >>> :{ + let + p :: M + p = fromList 5 5 + [ 0, 1, 0, 0, 0 + , 0, 0, 0.5, 0.5, 0 + , 0, 0, 0.5, 0.25, 0.25 + , 0, 0.25, 0, 0, 0.75 + , 1, 0, 0, 0, 0 + ] +:} +-} + +------------------------------------------------------------------------ +-- * Number of occurrences of a state in a test case + +-- | The fundamental matrix for absorbing chains. Its (i, j)-th entry is +-- the expected number of occurrences of state j prior to absorption at +-- the sink, given that one starts in state i. So the first row +-- indicates the expected occurence of each state starting from the +-- start state. +-- +-- >>> fundamental (minorMatrix 5 5 p :: M) +-- ┌ ┐ +-- │ 1.0 1.2307692307692306 1.2307692307692308 0.9230769230769231 │ +-- │ 0.0 1.2307692307692306 1.2307692307692308 0.9230769230769231 │ +-- │ 0.0 0.15384615384615385 2.1538461538461537 0.6153846153846154 │ +-- │ 0.0 0.3076923076923077 0.3076923076923077 1.2307692307692308 │ +-- └ ┘ +fundamental + :: M -- ^ Reduced matrix (n x n). + -> M -- ^ n x n +fundamental m = case inverse (identity (nrows m) .- m) of + Left err -> error $ "fundamental: " ++ err + Right r -> r + +-- | Expected variance of the occurrence for each state. The (i, j)-th +-- entry should be read in the same away as that of the fundamental +-- matrix above. +-- +-- >>> occVariance (fundamental (minorMatrix 5 5 p :: M)) +-- ┌ ┐ +-- │ 0.0 0.28402366863905315 2.5562130177514795 0.4970414201183433 │ +-- │ 0.0 0.28402366863905315 2.5562130177514795 0.4970414201183433 │ +-- │ 0.0 0.20118343195266267 2.4852071005917153 0.5207100591715977 │ +-- │ 0.0 0.3550295857988165 0.9230769230769231 0.2840236686390534 │ +-- └ ┘ +occVariance + :: M -- ^ Reduced matrix (n x n). + -> M -- n x n +occVariance n = n * (2 ^* diagonal 0 (getDiag n) - identity dim) - (n .* n) + where + dim = nrows n + +------------------------------------------------------------------------ +-- * Computing the long-run state probabilities + +-- | The Perron eigenvector (long-run occupancies/probabilities of states). +-- +-- >>> perron p +-- [0.1857142857142857,0.22857142857142854,0.22857142857142856,0.17142857142857143,0.1857142857142857] +perron :: M -- ^ Transition matrix (n x n). + -> V +perron p = Vector.map (/ l) n1 + where + dim = nrows p + + n1 :: V + n1 = getRow 1 (fundamental (minorMatrix dim dim p)) `Vector.snoc` 1 + + l :: Double + l = getElem 1 1 $ rowVector n1 * colVector (getDiag (identity dim)) + +------------------------------------------------------------------------ +-- * Sensitivity analysis + +-- XXX: Algorithm on p. 31 in report. +_sensitivities :: M -> M +_sensitivities = undefined + +------------------------------------------------------------------------ +-- * Other long run statistics + +-- | Compute the stimulus long-run occupancy. + +{- | >>> :{ +let s :: M + s = fromList 4 5 + [ 1, 0, 0, 0, 0 + , 0, 0.5, 0.5, 0, 0 + , 0, 0.5, 0.25, 0.25, 0 + , 0.25, 0, 0, 0.5, 0.25 + ] +in sigma s (perron p) +:} +[0.2807017543859649,0.2807017543859649,0.21052631578947367,0.17543859649122806,5.2631578947368425e-2] +-} +sigma :: M -- ^ Stimulus matrix (n-1 x n). + -> V -- ^ Perron eigenvector. + -> V +sigma stimulus pi = Vector.fromList + [ 1 / (1 - pi Vector.! (pred n)) + * sum [ pi Vector.! (pred i) * getElem i k stimulus + | i <- [1 .. (pred n)] + ] + | k <- [1..n] + ] + where + n = Vector.length pi + +-- ------------------------------------------------------------------------ +-- -- * Probability of occurrence for states + +-- | Compute the probability of occurrence for each state. +-- +-- >>> getRow 1 (nodeProbabilities p) +-- [1.0,1.0,0.5714285714285715,0.75] +nodeProbabilities + :: M -- ^ Transition matrix (n x n). + -> M -- n-1 x n-1 +nodeProbabilities p = n * (diagonal 0 (getDiag (1 ^/ n))) + where + n = fundamental (minorMatrix dim dim p) + dim = nrows p + +-- (Algorithm 9 on p. 34.) + +-- ------------------------------------------------------------------------ + +-- | Expected test case length. +-- +-- >>> expectedLength (fundamental (minorMatrix 5 5 p :: M)) +-- 4.384615384615385 +expectedLength + :: M -- ^ Fundamental matrix (n x n). + -> Double +expectedLength = norm_1 . rowVector . (getRow 1) + +-- ------------------------------------------------------------------------ + +-- | Success rate matrix. +-- +-- >>> expectedArcReliability Nothing (fromList 2 2 [10,10,10,10], fromList 2 2 [0,1,2,3]) +-- (┌ ┐ +-- │ 0.9166666666666666 0.8461538461538461 │ +-- │ 0.7857142857142857 0.7333333333333333 │ +-- └ ┘,┌ ┐ +-- │ 5.876068376068376e-3 9.298393913778529e-3 │ +-- │ 1.1224489795918367e-2 1.2222222222222223e-2 │ +-- └ ┘) +-- >>> :{ +-- expectedArcReliability +-- (Just (fromList 2 2 [10,10,10,10], fromList 2 2 [1,1,1,1])) +-- (fromList 2 2 [10,10,10,10], fromList 2 2 [0,1,2,3]) +-- :} +-- (┌ ┐ +-- │ 0.9130434782608695 0.875 │ +-- │ 0.84 0.8076923076923077 │ +-- └ ┘,┌ ┐ +-- │ 3.3081285444234404e-3 4.375e-3 │ +-- │ 5.169230769230769e-3 5.752794214332676e-3 │ +-- └ ┘) +expectedArcReliability + :: Maybe (M, M) + -> (M, M) + -> (M, M) -- ^ (mean, variance) +expectedArcReliability mprior (obsSuccs, obsFails) = + ( alpha ./ (alpha + beta) + , (alpha .* beta) ./ (ab .* ab .* (ab + ones)) + ) + where + m = nrows obsSuccs + n = ncols obsSuccs + + ones :: M + ones = Matrix.fromList m n [1,1..] + + -- In case no prior information is available, a_{i,j} = b_{i,j} = 1. + priorSuccs, priorFails :: M + (priorSuccs, priorFails) = + maybe (ones, ones) (bimap (+ ones) (+ ones)) mprior + + alpha, beta :: M + alpha = priorSuccs + obsSuccs + beta = priorFails + obsFails + + ab :: M + ab = alpha + beta + +transientReliability + :: M -- ^ Reduced transition matrix (n-1 x n). + -> Maybe (M, M) -- ^ Prior successes and failures (n-1 x n). + -> (M, M) -- ^ Observed successes and failures (n-1 x n). + -> (M, M) -- ^ Reliability vector and reliability variance vector. +transientReliability q mprior obs = + (rstar, vstar) + where + n = ncols q + + (mean, var) = expectedArcReliability mprior obs + + r1 :: M + r1 = mean + + fancyR :: M + fancyR = q .* r1 + + fancyRdot :: M -- n-1 x n-1 + fancyRdot = submatrix 1 (pred n) 1 (pred n) fancyR + + w :: M -- n-1 x 1 + w = colVector $ getCol n fancyR + + rstar :: M -- n-1 x 1 + rstar = fundamental fancyRdot * w + + r2 :: M + r2 = r1 .* r1 + var + + fancyR' :: M -- n-1 x n-1 + fancyR' = q .* r2 + + fancyRdot' :: M -- n-1 x n-1 + fancyRdot' = submatrix 1 (pred n) 1 (pred n) fancyR' + + w' :: M -- n-1 x 1 + w' = colVector $ getCol n fancyR' + + rstar' :: M -- n-1 x 1 + rstar' = fundamental fancyRdot' * w' + + vstar = rstar' - rstar .* rstar + +singleUseReliability + :: M -- ^ Reduced transition matrix (n-1 x n). + -> Maybe (M, M) -- ^ Prior successes and failures (n-1 x n). + -> (M, M) -- ^ Observed successes and failures (n-1 x n). + -> (Double, Double) +singleUseReliability q mprior obs = + bimap (head . toList) (head . toList) $ transientReliability q mprior obs + +------------------------------------------------------------------------ + +{- +singleUseReliabilityIO :: (KnownNat n, KnownNat (n - 1)) + => (1 <= n - 1, n - 1 <= n, (n - (n - 1)) ~ 1) + + => L (n - 1) n -- ^ Transition matrix without transitions to + -- the sink. + + -> FilePath -- ^ Filepath where prior successes are/will be + -- saved. + + -> FilePath -- ^ Filepath where prior failures are/will be + -- saved. + + -> (L (n - 1) n, L (n - 1) n) -- ^ Observed successes and failures. + -> IO Double +singleUseReliabilityIO q fps fpf (s, f) = do + mprior <- load2StaticMatrices fps fpf + case mprior of + Nothing -> do + -- We need the elementwise @max 1@ below because otherwise we get division + -- by zero in 'expectedArcReliability'. For details see the 2017 paper by by L. Lin, + -- Y. Xue and F. Song in the README. + saveStaticMatrix fps "%.1f" (dmmap (max 1) s) + saveStaticMatrix fpf "%.1f" (dmmap (max 1) f) + Just (ps, pf) -> do + saveStaticMatrix fps "%.1f" (ps + s) + saveStaticMatrix fpf "%.1f" (pf + f) + return (singleUseReliability q mprior (s, f)) + +loadStaticMatrix :: (KnownNat m, KnownNat n) => FilePath -> IO (Maybe (L m n)) +loadStaticMatrix = fmap (join . fmap create) . D.loadMatrix' + +load2StaticMatrices :: (KnownNat m, KnownNat n, KnownNat o, KnownNat p) + => FilePath -> FilePath -> IO (Maybe (L m n, L o p)) +load2StaticMatrices fp1 fp2 = liftMA2 (loadStaticMatrix fp1) (loadStaticMatrix fp2) + where + liftMA2 :: (Monad m, Applicative f) => m (f a) -> m (f b) -> m (f (a, b)) + liftMA2 = liftM2 (liftA2 (,)) + +saveStaticMatrix :: (KnownNat m, KnownNat n) + => FilePath + -> String -- ^ "printf" format (e.g. "%.2f", "%g", etc.) + -> L m n -> IO () +saveStaticMatrix fp format = D.saveMatrix fp format . unwrap +-} + +------------------------------------------------------------------------ + +-- | Kullback-Leibler matrix discrimination. + +{- | >>> :{ +let + s = fromList 4 5 + [ 1, 0, 0, 0, 0 + , 0, 0.5, 0.5, 0, 0 + , 0, 0.5, 0.25, 0.25, 0 + , 0.25, 0, 0, 0.5, 0.25 + ] + es = Vector.fromList [4, 5, 7, 3, 4] + et = fromList 4 5 + [ 4, 0, 0, 0, 0 + , 0, 3, 2, 0, 0 + , 0, 4, 1, 2, 0 + , 1, 0, 0, 1, 1 + ] +in kullbackLeibler p s es et +:} +3.440540391434413e-2 +-} +kullbackLeibler + :: M -- ^ Transitions (n x n). + -> M -- ^ Stimulis (m x n). + -> V -- ^ State visitations (n). + -> M -- ^ Stimulus execution (m x n). + -> Double -- ^ Discrimination. +kullbackLeibler p s es et = Vector.sum k' + where + m = nrows s + n = ncols p + + pi' :: V + pi' = Vector.take m (perron p) + + es' :: V + es' = Vector.take m es + + es'' :: M + es'' = let v = colVector es' + in List.foldl' (<|>) v (take (pred n) (repeat v)) + + t :: M + t = es'' ./ et + + k :: V + k = Vector.fromList $ map sum $ toLists $ (1 / log 2) ^* (s .* fmap log' (s .* t)) + + k' :: V + k' = Vector.zipWith (*) pi' k + + log' :: Double -> Double + log' x + | x == 0 = 0 + | isNaN x = 0 + | otherwise = log x diff --git a/test/Spec.hs b/test/Spec.hs new file mode 100644 --- /dev/null +++ b/test/Spec.hs @@ -0,0 +1,1 @@ +{-# OPTIONS_GHC -F -pgmF tasty-discover -optF --tree-display #-} diff --git a/test/Unit.hs b/test/Unit.hs new file mode 100644 --- /dev/null +++ b/test/Unit.hs @@ -0,0 +1,118 @@ +module Unit + ( unit_expectExpectedArcReliability + , unit_expectTransientReliability + , unit_occurenceMean + , unit_occurenceVar + , unit_docTest + ) where + +import Data.Matrix +import qualified Data.Vector as Vector +import Prelude hiding + (pi) +import Test.DocTest + (doctest) +import Test.Tasty.HUnit + (Assertion, (@?), (@?=)) + +import MarkovChain + +------------------------------------------------------------------------ + +-- cf. A Simpler and More Direct Derivation of System Reliability Using Markov Chain Usage Models (2017) +-- by Lan Lin, Yufeng Xue and Fengguang Song + +-- Transient usage model (transitions from sink omitted). +q :: M +q = fromList 2 3 + [ 0, 0.5, 0.5 + , 0, 0.5, 0.5 + ] + +successes :: M +successes = fromList 2 3 + [ 1, 2, 3 + , 4, 5, 6 + ] + +failures :: M +failures = fromList 2 3 + [ 1, 0, 1 + , 0, 1, 0 + ] + +-- Observed transient success rate. +sr :: M +(sr, _) = expectedArcReliability Nothing (successes, failures) + +unit_expectExpectedArcReliability :: Assertion +unit_expectExpectedArcReliability = sr @?= expected + where + expected :: M + expected = fromList 2 3 + [ 2/4, 3/4, 4/6 + , 5/6, 6/8, 7/8 + ] + +-- Transient reliability matrix. +tr :: M +tr = fst $ transientReliability q Nothing (successes, failures) + +unit_expectTransientReliability :: Assertion +unit_expectTransientReliability = + (norm_F (tr - expected)) <= 1.0e-6 @? "differs from expected" + where + a = (4/6)/2 + b = (3/4)*(7/8)/4 + c = (7/8)/2 + d = 1 - (6/8)/2 + + expected :: M + expected = fromList 2 1 + [ a + b/d + , c/d + ] + +------------------------------------------------------------------------ + +-- cf. Computations for Markov Chain Usage Models (2000) +-- by S. J. Prowell + +-- Usage model. +p :: M +p = fromList 5 5 + [ 0, 1, 0, 0, 0 + , 0, 0, 0.5, 0.5, 0 + , 0, 0, 0.5, 0.25, 0.25 + , 0, 0.25, 0, 0, 0.75 + , 1, 0, 0, 0, 0 + ] + +q' :: M +q' = minorMatrix 5 5 p + +-- State occurences. +unit_occurenceMean :: Assertion +unit_occurenceMean = + norm_F (rowVector actual .- rowVector expected) <= 1.0e-3 @? "differs from expected" + where + actual :: V + actual = getRow 1 (fundamental q') + + expected :: V + expected = Vector.fromList + [ 1.0, 1.231, 1.231, 0.9231 ] + +unit_occurenceVar :: Assertion +unit_occurenceVar = + norm_F (rowVector actual .- rowVector expected) <= 1.0e-3 @? "differs from expected" + where + actual :: V + actual = getRow 1 (occVariance (fundamental q')) + + expected :: V + expected = Vector.fromList + [ 0, 0.284, 2.556, 0.497 ] + +unit_docTest :: IO () +unit_docTest = doctest ["src/MarkovChain.hs"]