1/*
    2Hidden Markov model for modeling DNA sequences.
    3The model has three states, q1, q2 and end, and four output symbols,
    4a, c, g, and t, corresponding to the four nucleotides (letter).
    5From
    6Christiansen, H. and Gallagher, J. P. 2009. Non-discriminating arguments and
    7their uses. In International Conference on Logic Programming. LNCS,
    8vol. 5649. Springer, 55-69.
    9
   10*/
   11:- use_module(library(pita)).   12
   13:- if(current_predicate(use_rendering/1)).   14:- use_rendering(c3).   15:- use_rendering(graphviz).   16:- use_rendering(table,[header(['Multivalued variable index','Rule index','Grounding substitution'])]).   17:- endif.   18
   19:- pita.   20
   21:- begin_lpad.   22
   23% hmm(O): O is the output sequence
   24% hmm1(S,O): O is the output sequence and S is the sequence of states
   25% hmm(Q,S0,S,O):  from state Q and previous state S0, generates output O and
   26% sequence of states S
   27
   28hmm(O):-hmm1(_,O).
   29% O is an output sequence if there is a state seuqnece S such that hmm1(S,O)
   30% holds
   31
   32hmm1(S,O):-hmm(q1,[],S,O).
   33% O is an output sequence and S a state sequence if the chain stats at state
   34% q1 and ends generating state sequence S and output sequence O
   35
   36hmm(end,S,S,[]).
   37% an HMM in state end terminates the sequence without emitting any symbol
   38
   39hmm(Q,S0,S,[L|O]):-
   40	Q\= end,
   41	next_state(Q,Q1,S0),
   42	letter(Q,L,S0),
   43	hmm(Q1,[Q|S0],S,O).
   44% an HMM in state Q different from end goes in state Q1, emits the letter L
   45% and continues the chain
   46
   47map_query next_state(q1,q1,S):1/2;next_state(q1,q2,S):0.45;next_state(q1,end,S):0.05.
   48% from state q1 the HMM can go to q1, q2 or end with equal probability
   49
   50map_query next_state(q2,q1,S):0.45;next_state(q2,q2,S):1/2;next_state(q2,end,S):0.05.
   51% from state q2 the HMM can go to q1, q2 or end with equal probability
   52
   53
   54map_query letter(q1,a,S):0.4;letter(q1,c,S):0.3;letter(q1,g,S):0.2;letter(q1,t,S):0.1.
   55% from state q1 the HMM emits one of the letters with equal probability
   56
   57map_query letter(q2,a,S):0.1;letter(q2,c,S):0.2;letter(q2,g,S):0.3;letter(q2,t,S):0.4.
   58% from state q1 the HMM emits one of the letters with equal probability
   59
   60:- end_lpad.   61
   62
   63state_diagram(digraph(G)):-
   64    findall(edge((A -> B),[label=P]),
   65      (clause(next_state(A,B,_,_,_),
   66        (get_var_n(_,_,_,_,Probs,_),equality(_,_,N,_))),
   67        nth0(N,Probs,P)),
   68      G).