Computing the high-end distribution from lazy permutations

Searching matrix paths using the built-in laziness of Haskell

warning: source code ahead!

Finding all paths through a matrix

Calculate all selections of one element from each row

..by combining each element from the first row with all paths through the rest of the matrix

paths :: [[a]] -> [[a]]

paths [r] = map return r

paths (r:rs) = [x:y | y <- paths rs, x <- r]

Finding paths in descending order

As before, generate the list of combinations of

previous lists

elements from the current row

but merge-sort them (by the sum of elements)

paths [r] = map return $ reverse $ sort r

paths (r:rs) = foldr1 merge

[[x:y | y <- paths rs] | x <- r]

From paths to scored probabilites

Not really interested in the paths (strings)

Must add probabilities for each character

Some type definitions:

type PSSM = [[Score]]

type PVector = [Prob]

type Score = Int

type Prob = Double

From paths to scored probabilites

A PVector is a list of probabilities for the given alphabet

Same order (and length!) as the rows

prep_row :: PVector -> [Score] -> [(Score,Prob)]

prep_row p r = reverse $ sort $ zip r p

From paths to scored probabilities

Only retain score and prob for each path

Rewrite merge to work on score-prob pairs instead of lists of scores

paths p [r] = prep_row p r

paths p (r:rs) = foldr1 merge

[[combine x y | y <- paths p rs] | x <- prep_row p r]

combine (s1,p1) (s2,p2) = (s1+s2,p1*p2)

Applying Dynamic Programming

Paths with the same score can be combined

Merge-sorting brings similar scores together

redund ((s,p):(s',p'):xs) =

if s /= s' then (s,p) : redund ((s',p'):xs)

else redund ((s,p+p'):xs)

redund xs = xs

merge' xs = redund . merge xs

Timing 1000 PRINTS matrices

Complete set of paths:

Approx 80 min.

PoSSuM ~ 10x faster

Specific p-values:

1e-10 38s (1.2x)

1e-20 104s (2.5x)

1e-30 369s (6x)

Conclusions

Haskell's laziness lets us implent the probability calculation in a simple and straightforward manner

Lags speedwise behind the C implementation by an order of magnitude

2004-06-15T23:56:16