# Rosalind: BA10A — Compute the Probability of a Hidden Path
# https://rosalind.info/problems/ba10a/
#
# Given: A hidden path pi, the states of an HMM, and its transition matrix.
# Return: Pr(pi).

let hidden_path = "AABBBAABABAAAABBBBAABBABABBBAABBAAAABABAABBABABBAB"

# No emissions are involved, so the model needs only its states and how they
# follow one another.
let model = {
    states: ["A", "B"],
    transition: {
        A: { A: 0.194, B: 0.806 },
        B: { A: 0.273, B: 0.727 },
    },
}

# Every state is equally likely to start, then each step multiplies in one
# transition. Fifty of them take the answer down to 1e-19, which is what makes
# the log-space treatment in BA10C and BA10D necessary rather than fussy.
let probability = hmm_path_probability(hidden_path, model)

# Printed as a mantissa and an exponent; the plain decimal expansion of 1e-19
# is unreadable next to the published answer.
println("Result:   " + str(round(probability * 1e19, 6)) + "e-19")
println("Expected: 5.017329e-19")

fn test_ba10a_hidden_path_probability() {
    assert abs(probability - 5.01732865318e-19) < 1e-30,
        "BA10A: got " + str(probability)
    # The same thing computed by hand, to check the builtin agrees with the
    # definition rather than only with itself.
    let steps = chars(hidden_path)
    let by_hand = range(1, len(steps))
        |> reduce(|running, i| running * model.transition[steps[i - 1]][steps[i]], 0.5)
    assert abs(probability - by_hand) < 1e-30, "BA10A: builtin and hand computation disagree"
}