F#. Lazy Levenshtein Distance.

Update : Hier die Visualisierung.

Jetzt ist mir endlich gelungen den Lazy Levenshtein Distance Algorithmus aus dem letzten Post in F# zu implementieren.

// levenshtein.fs
// According to the article, the worst-case complexity is O(|A|*|B|).
// In this case, the array version is much faster, which is probably due to the overhead 
// of the use of Lazy and LazyList delayed methods in F#.
// However, if A=B the complexity is now O(|A|) because only the main diagonal is evaluated.
// In this case, lazy version is faster.
    //from http://www.haskell.org/haskellwiki/Edit_distance.
    // "An entry depends on three neighbours which lie on the diagonal below, the current diagonal and the diagonal above.
    //  Each diagonal therefore depends on the diagonal below and the diagonal above where a row depends only on the row above"    
    let inline editDist sa sb =
        let min3 x y z = if x < y then x else min y (LazyList.head z)
        let lab = List.length sa - List.length sb
        let rec mainDiag:_ Lazy = 
            lazy(oneDiag sa sb (LazyList.delayed (fun () -> LazyList.head uppers.Value)) 
                               (LazyList.consDelayed -1 (fun ()-> LazyList.head lowers.Value)))
        //upper diagonals
        and uppers : _ Lazy = lazy(eachDiag sa sb (LazyList.consDelayed (mainDiag.Value) (fun ()-> uppers.Value)))
        //lower diagonals. note swap sb sa !
        and lowers : _ Lazy = lazy(eachDiag sb sa (LazyList.consDelayed (mainDiag.Value) (fun ()-> lowers.Value)))
        // 'a list -> 'a list -> LazyList<LazyList<int>>
        and eachDiag a b diag = 
            match a, b, diag with
            | _, [], _ -> LazyList.empty
            | a, (bch :: bs), ( LazyList.Cons (lastDiag, diags)) -> 
                let nextDiag = LazyList.delayed (fun () -> LazyList.head (LazyList.tail diags))
                LazyList.consDelayed (oneDiag a bs nextDiag lastDiag) (fun ()-> (eachDiag a bs diags))
        // 'a list -> 'a list -> LazyList<int> -> LazyList<int> -> LazyList<int>
        and oneDiag a b diagAbove diagBelow  = 
            // nw - north-west, n - north, w - west.
            // 'a list -> 'a list -> int -> LazyList<int> -> LazyList<int> -> LazyList<int>
            let rec doDiag a b nw n w = 
                match a, b with
                | [], _ -> LazyList.empty
                | _, [] -> LazyList.empty 
                | (ach :: achs), (bch :: bchs) -> 
                    
                    let me  = if ach = bch then nw else 1 + min3 (LazyList.head w) nw n
                    LazyList.consDelayed me (fun () -> 
                                                doDiag achs bchs me                         // hope these
                                                    (LazyList.delayed (fun () -> LazyList.tail n))  // <---    
                                                    (LazyList.delayed (fun () -> LazyList.tail w))) // <--- not evaluated.
               
            let firstelt = 1 + (LazyList.head diagBelow)
            LazyList.consDelayed firstelt (fun ()-> doDiag a b firstelt diagAbove (LazyList.tail diagBelow))

        if lab = 0      then mainDiag.Value
        else if lab > 0 then (LazyList.toArray lowers.Value).[lab - 1]
        else                 (LazyList.toArray (uppers.Value)).[-1 - lab]

    let inline lazyDist (s1:string) (s2:string) =
        match s1.Length,s2.Length with
        |0, l2 -> l2
        |l1, 0 -> l1
        | _ -> 
            let res = editDist (List.ofSeq s1) (List.ofSeq s2) |> LazyList.toArray
            res.[res.Length - 1]

let inline timeExec f a b s =
    let timer = new System.Diagnostics.Stopwatch()
    timer.Start()
    let res = f a b
    timer.Stop()
    printfn "%A." s
    printfn "distance = %A: Ellapsed Time: %A ticks, %A ms." res timer.ElapsedTicks timer.ElapsedMilliseconds

let str1 = String.replicate 500 "abcd" 
let str2 = String.replicate 500 "defg"  
let str3 = String.replicate 1000 "a"  
let str4 = (String.replicate 20 "aba")+(String.replicate 500 "aa") + "aaa" + (String.replicate 500 "aa") + (String.replicate 20 "aba") 
let str5 = (String.replicate 20 "aca")+(String.replicate 500 "aa") + "bbb" + (String.replicate 500 "aa") + (String.replicate 20 "aca")

[0..10] |> List.map (fun _ -> 
    printfn "-----------------------"
    timeExec levenshteinDistanceArray str1 str2 "levenshtein Distance with Array. worst-case."
    timeExec lazyDist str1 str2 "lazy levenshtein Distance. worst-case."
    timeExec levenshteinDistanceArray str3 str3 "levenshtein Distance with Array. special case of similar strings."
    timeExec lazyDist str3 str3 "lazy levenshtein Distance. special case of similar strings."
    timeExec levenshteinDistanceArray str4 str5 "levenshtein Distance with Array."
    timeExec lazyDist str4 str5 "lazy levenshtein Distance."
    )|>ignore

-----------------------
"levenshtein Distance with Array. worst-case.".
distance = 1502: Ellapsed Time: 3858477L ticks, 269L ms,
"lazy levenshtein Distance. worst-case.".
distance = 1502: Ellapsed Time: 147370307L ticks, 10292L ms.
"levenshtein Distance with Array. special case of similar strings."
distance = 0: Ellapsed Time: 1072085L ticks, 74L ms.
"lazy levenshtein Distance. special case of similar strings.".
distance = 0: Ellapsed Time: 115395L ticks, 8L ms.
"levenshtein Distance with Array.".
distance = 43: Ellapsed Time: 4217920L ticks, 294L ms.
"lazy levenshtein Distance.".
distance = 43: Ellapsed Time: 3744706L ticks, 261L ms.
-----------------------
"levenshtein Distance with Array. worst-case.".
distance = 1502: Ellapsed Time: 3896606L ticks, 272L ms.
"lazy levenshtein Distance. worst-case.".
distance = 1502: Ellapsed Time: 146242131L ticks, 10213L ms.
"levenshtein Distance with Array. special case of similar strings."
distance = 0: Ellapsed Time: 1038478L ticks, 72L ms.
"lazy levenshtein Distance. special case of similar strings.".
distance = 0: Ellapsed Time: 10024L ticks, 0L ms.
"levenshtein Distance with Array.".
distance = 43: Ellapsed Time: 4204501L ticks, 293L ms.
"lazy levenshtein Distance.".
distance = 43: Ellapsed Time: 3610026L ticks, 252L ms.
-----------------------
"levenshtein Distance with Array. worst-case.".
distance = 1502: Ellapsed Time: 4039147L ticks, 282L ms.
"lazy levenshtein Distance. worst-case.".
distance = 1502: Ellapsed Time: 145049786L ticks, 10130L ms.
"levenshtein Distance with Array. special case of similar strings."
distance = 0: Ellapsed Time: 1180816L ticks, 82L ms.
"lazy levenshtein Distance. special case of similar strings.".
distance = 0: Ellapsed Time: 27104L ticks, 1L ms.
"levenshtein Distance with Array.".
distance = 43: Ellapsed Time: 4247992L ticks, 296L ms.
"lazy levenshtein Distance.".
distance = 43: Ellapsed Time: 3652880L ticks, 255L ms.

Update
Write "min3" and "let me = if ach = bch then nw else 1 + min3 (LazyList.head w) nw n" explicitly
in order to avoid the unnecessary delay.
give an additional performance gain.

// levenshtein.fs 
    let inline editDist sa sb =
        //let min3 x y z = if x < y then x else min y (LazyList.head z)
        let lab = List.length sa - List.length sb
        let rec mainDiag:_ Lazy = 
            lazy(oneDiag sa sb (LazyList.delayed (fun () -> LazyList.head uppers.Value)) 
                               (LazyList.consDelayed -1 (fun ()-> LazyList.head lowers.Value)))
        //upper diagonals
        and uppers : _ Lazy = lazy(eachDiag sa sb (LazyList.consDelayed (mainDiag.Value) (fun ()-> uppers.Value)))
        //lower diagonals. note swap sb sa !
        and lowers : _ Lazy = lazy(eachDiag sb sa (LazyList.consDelayed (mainDiag.Value) (fun ()-> lowers.Value)))
        // 'a list -> 'a list -> LazyList<LazyList<int>>
        and eachDiag a b diag = 
            match a, b, diag with
            | _, [], _ -> LazyList.empty
            | a, (bch :: bs), ( LazyList.Cons (lastDiag, diags)) -> 
                let nextDiag = LazyList.delayed (fun () -> LazyList.head (LazyList.tail diags))
                LazyList.consDelayed (oneDiag a bs nextDiag lastDiag) (fun ()-> (eachDiag a bs diags))
        // 'a list -> 'a list -> LazyList<int> -> LazyList<int> -> LazyList<int>
        and oneDiag a b diagAbove diagBelow  = 
            // nw - north-west, n - north, w - west.
            // 'a list -> 'a list -> int -> LazyList<int> -> LazyList<int> -> LazyList<int>
            let rec doDiag a b nw n w = 
                match a, b with
                | [], _ -> LazyList.empty
                | _, [] -> LazyList.empty 
                | (ach :: achs), (bch :: bchs) -> 
                    if (ach  = bch) then 
                //case if ach = bch then nw
                        LazyList.consDelayed nw (fun ()-> 
                                                    doDiag achs bchs nw 
                                                        (LazyList.delayed (fun ()-> LazyList.tail n)) 
                                                        (LazyList.delayed (fun ()-> LazyList.tail w)))
                // case else 1 + min3 (LazyList.head w) nw n
                    else if (LazyList.head w) < nw then
                    // case let min3 x y z = if x < y then x ...
                        let me = 1 + (LazyList.head w)
                        LazyList.consDelayed me (fun ()-> 
                                                    doDiag achs bchs me 
                                                        (LazyList.delayed (fun ()-> LazyList.tail n)) 
                                                        (LazyList.tail w))
                    // case let min3 x y z = ... else min y (LazyList.head z)
                    else
                        let me = 1 + min nw (LazyList.head n)
                        LazyList.consDelayed me (fun ()-> 
                                                    doDiag achs bchs me 
                                                        (LazyList.tail n)
                                                        (LazyList.tail w))               
            let firstelt = 1 + (LazyList.head diagBelow)
            LazyList.consDelayed firstelt (fun ()-> doDiag a b firstelt diagAbove (LazyList.tail diagBelow))

        if lab = 0      then mainDiag.Value
        else if lab > 0 then (LazyList.toArray lowers.Value).[lab - 1]
        else                 (LazyList.toArray (uppers.Value)).[-1 - lab]

    let inline lazyDist (s1:string) (s2:string) =
        match s1.Length,s2.Length with
        |0, l2 -> l2
        |l1, 0 -> l1
        | _ -> 
            let res = editDist (List.ofSeq s1) (List.ofSeq s2) |> LazyList.toArray
            res.[res.Length - 1]

F# Type-directed memoization. IntTrie, knapsack problem and levenshtein distance.

Nach wochenlangem Bewerbungsstress und Frustration darüber, dass es wohl kaum F#-Jobstellen auf dem Markt gibt und die wenigen, die da sind, richten sich ausschlislich an Hochschulabsolventen, habe ich mich weiter mit der "Type-directed Memoization" beschäftigt. Dabei endeckte ich eine interessante OCaml Implementation.
Ich suchte nach konkreten praktischen Beispielen und merkte sehr schnell, dass in diesem Zusammenhang oft zwei Algorithmen genannt werden. Das Knapsack Problem
und die Levenshtein-Distanz. Hier ein paar Links zu diesen Themen.

Solving the 0-1 knapsack problem using continuation-passing style with memoization in F#.
    Haskell Version.
    Haskell Version mit der IntTrie-Datenstruktur.
    Die Levenshtein-Distanz auf Rosetta Code Seite.
    "Haskell function computes the edit distance in O(length a * (1 + dist a b)) time complexity".
    Lazy Levenshtein Distanz.

Von IntTrie war ich sofort begeistert und versuchte die abgespeckte Version - also nur positive Integers - nach F# zu übertragen.

// TypeMemo.fs
namespace TypeMemoization
// from http://hackage.haskell.org/packages/archive/data-inttrie/0.0.7/doc/html/src/Data-IntTrie.html
module BitTrie =
    
    type BitTrie<'a> = BitTrie of Lazy<'a> * Lazy<BitTrie<'a>> * Lazy<BitTrie<'a>>
    // A trie from positiv integers to values of type a. 
    type IntTrie<'a> = IntTrie of Lazy<'a> * BitTrie<'a>
    
    let inline testBit x = (x &&& 1) <> 0

    let rec fmap f (BitTrie (x, l, r)) = 
        BitTrie( lazy(f x.Value), 
                 lazy(fmap f l.Value), 
                 lazy(fmap f r.Value) )

    let identityPositive = 
        let rec go x = 
            BitTrie (x, 
                     lazy(fmap (fun n-> n <<< 1) (go x)), 
                     lazy(fmap (fun n -> (n <<< 1) ||| 1) (go x)))
        go (lazy(1))

    let inline fmapi f (IntTrie(z, pos)) = 
        IntTrie(lazy(f z.Value), fmap f pos)

    //The identity trie.
    let identity = IntTrie (lazy(0), identityPositive)

    let inline toTrie f  =  fmapi f identity

    let rec applyPositive (BitTrie (one, even, odd)) x =
        match x with
        | i when i = 1 ->       one
        | i when testBit i  ->  applyPositive odd.Value (x >>> 1) 
        | otherwise   ->        applyPositive even.Value (x >>> 1)

    // Apply the trie to an argument.
    let inline apply (IntTrie(z, pos)) x =
        match x with
        | 0 -> z.Value
        |_ -> (applyPositive pos x).Value 
    
    let inline memo f = apply (toTrie f)

    // Memoize a two argument function (just apply the table directly for
    // single argument functions).
    let inline memo2 f = memo (memo << f)

//knapsack.fs
namespace TypeMemoization

module knapsack =
    open System

    let inline genItems n = 
        match n with
        | 0 -> Array.empty 
        | _ -> Array.init n 
                        ( fun i ->
                            let weight = i % 5
                            let value = (float)(weight * i)
                            weight, value )
    let inline max (x:float) (y:float) = max x y 

    let inline knapsackOriginal desiredWeight (items:_[]) =               
        let inline weightOf i = fst items.[i-1]
        let inline valueOf i = snd items.[i-1] 

        let rec knapsack' i w  = 
            match i, w  with
            | 0, _ | _, 0 -> 0.
            | i, w    -> 
                match i with
                | i' when (weightOf i') > w -> 
                    knapsack' (i' - 1) w           
                | _ -> 
                    max (knapsack' (i - 1) w)  ((knapsack' (i - 1) (w - weightOf i)) + valueOf i)
                    
        knapsack' items.Length desiredWeight

    let inline knapsack weight (value:int->float) =                
        
        let rec knapsack' i w  = 
            match i, w  with
            | 0, _ | _, 0 -> 0.
            | i, w    -> 
                match i with
                | i' when (weight i') > w ->
                    knapsackMemo (i' - 1) w           
                | _ -> 
                    max (knapsackMemo (i - 1) w)  ((knapsackMemo (i - 1) (w - weight i)) + value i)
                    
        and knapsackMemo  = BitTrie.memo2 knapsack'
        knapsackMemo  

    let inline knapsackMemoized desiredWeight (items:_[]) =
        let inline weightOf i = fst items.[i-1]
        let inline valueOf i = snd items.[i-1] 

        knapsack weightOf valueOf (items.Length) (desiredWeight)

Zwar ist die "memoized" Version schneller als der originale Knapsack-Algorithmus, aber leider viel langsamer als Zach Bray's Version und folglich um X-faches langsamer als die imperative Variante. Es kann aber auch sein, dass meine F# Implementierung von IntTrie nicht die beste ist.

Mit der Levenshtein Distanz sieht es noch schlimmer aus. Wieder ist die "memoized" Variante besser als die naive Version, aber mit der Array-Lösung kann sie überhaupt nicht mithalten. Interessant ist der bereits erwähnte Lazy Levenshtein Algorithmus. Dieser in F# zu übertragen ist mir leider nicht gelungen.
Also ist mein Frust nur noch tiefer geworden.

// from http://research.microsoft.com/en-us/um/people/simonpj/papers/assoc-types/fun-with-type-funs/typefun.pdf.
// chapter "3.1 Type-directed memoization".
namespace TypeMemoization
module BitTrie =
    ...
module TypeMemo =
    open BitTrie

    type IFromTable<'a,'w > =
        abstract inline fromTable : 'a->'w
    
    let inline fromTable t = (t :> IFromTable<_,_>).fromTable
    // "we can memoise any function from Bool by storing its two
    // return values as a lazy pair. This lazy pair is the memo table."
    type BoolTable<'w> = 
        | BTable of Lazy<'w> * Lazy<'w> 
        interface IFromTable<bool,'w> with
            member inline x.fromTable b = 
                match x with
                | BTable(x,y) -> if b then x.Force() else y.Force()
    
    let inline boolToTable f = 
        BTable (lazy(f true), lazy(f false))
    
    
    //"memoise functions from any sum type, such as the type Either."
    type Either<'a,'b>= 
        |Left of 'a
        |Right of 'b
    // "We can memoise a function from Either a b by storing a lazy pair of a
    // memo table from a and a memo table from b. That is, we take advantage
    // of the isomorphism between the function type Either a b -> w and the
    // product type (a -> w, b -> w)."
    type DiscriminatedUnionTable<'a,'b,'w> = 
        | STable of IFromTable<'a,'w> * IFromTable<'b,'w> 
        interface IFromTable<Either<'a,'b>,'w> with
            member inline x.fromTable e = 
                match x, e with
                | STable (t,_), Left  v   -> t.fromTable v
                | STable (_,t), Right v   -> t.fromTable v
    
    let inline discriminatedUnionToTable f fa fb =
        STable (fa (f<<Left),fb (f<<Right))
    
    
    // "Dually, we can memoise functions from the product type (a,b) by storing a memo table
    // from a whose entries are memo tables from b."
    type ProductTable<'a,'b,'w> = 
        | PTable of IFromTable<'a, IFromTable<'b,'w>>
        interface IFromTable<('a*'b),'w> with
            member inline x.fromTable p = 
                match x, p with
                | PTable t,(a,b)-> (t.fromTable a).fromTable b

    let inline productToTable f fa fb = 
        PTable ((fa (fun a -> fb (fun b -> f (a, b)) :> IFromTable<_,_> )) :> IFromTable<_,_>)   
   

    // "A list is a combination of a sum, a product, and recursion.
    // Since a list is either empty or not, ListTable<'a,'w> is represented by a pair, whose first component is the result of applying
    // the memoised function f to the empty list, and whose second component
    // memoises applying f to non-empty lists."
    type ListTable<'a,'w> =
        | LTable of Lazy<'w> * IFromTable<'a, IFromTable<'a list, 'w>>
        interface IFromTable<'a list,'w> with
            member inline x.fromTable l =
                match x, l with
                | LTable(t, _), [] -> t.Force() 
                | LTable(_, t), x :: xs -> fromTable (fromTable t x ) xs
                    
    let rec listToTable f fa =
        LTable (lazy(f []), fa (fun x -> listToTable (fun xs ->  f (x::xs)) fa :> IFromTable<_,_>) :> IFromTable<_,_>)
    let inline flip1 f a b c = f c a b

    type CharTable<'w> =
        | CharTable of IntTrie<'w>  
        interface IFromTable<char, 'w> with
            member inline x.fromTable c =
                match x with
                | CharTable t ->  apply t (int c)
    let inline charToTable f = 
        CharTable (toTrie (f << char))
    
    let inline memoCharList2 f = 
        let memo g =
            listToTable g charToTable  |> fromTable
        memo (memo << f)

Levenshtein Distanz.

//levenshtein .fs
namespace TypeMemoisation
module levenshtein =
    open System
    
    let inline private naiveLevenshteinDistance del sub ins  =  
        let rec inner s1 s2 =
            match s1, s2 with
            | s1,     []     -> ins * List.length s1 
            | [],     s2     -> ins * List.length s2 
            | x :: xs, y :: ys ->
                match x = y with
                | true -> inner xs ys
                | _ -> List.min [ del + inner xs s2; 
                                 sub + inner s1 ys; 
                                 ins + inner xs ys] 
        inner

    let inline runNaiveLevenshteinDistance (s1:string) (s2:string) =
        naiveLevenshteinDistance 1 1 1 (List.ofSeq s1) (List.ofSeq s2)
    //memoized version.
    let inline private levenshteinDistance del sub ins  =  
        let rec inner s1 s2 =
            match s1, s2 with
            | s1,     []     -> ins * List.length s1 
            | [],     s2     -> ins * List.length s2 
            | x :: xs, y :: ys ->
                match x = y with
                | true -> memo xs ys
                | _ -> min (del + memo xs s2)   
                                 (min (sub + memo s1 ys) (ins + memo xs ys)) 
        and memo = 
                TypeMemo.memoCharList2 inner
        memo

    let inline levenshteinDistanceMemoized (s1:string) (s2:string) =
        levenshteinDistance 1 1 1 (List.ofSeq s1) (List.ofSeq s2)
    
    //array version.
    let inline levenshteinDistanceArray (s1:string) (s2:string) =
         let sa, sb:char [] * char [] = s1.ToCharArray(), s2.ToCharArray()
         let len = Array.length sa
         let m = len - 1

         let inline compute z xc = min (z+1) xc

         let inline transform (narr : int []) chb =
              Array.zip3 sa.[..m] narr.[..m] narr.[1..m+1] 
              |> Array.map (fun (cha, x, y) -> min (y + 1) (x + abs(compare cha chb)) )
              |> Array.scan compute (narr.[0] + 1)
         let result = Array.fold transform [|0..len|] sb
         result.[result.Length - 1]