wolffd@0: function bnet = fixed_args_mk_motif_hhmm(motif_length, motif_pattern, background_char)
wolffd@0: %
wolffd@0: % BNET = MK_MOTIF_HHMM(MOTIF_LENGTH)
wolffd@0: % Make the following HHMM
wolffd@0: %
wolffd@0: %    S2 <----------------------> S1
wolffd@0: %    |                           |
wolffd@0: %    |                           |
wolffd@0: %   M1 -> M2 -> M3 -> end        B1 -> end
wolffd@0: %
wolffd@0: % where Mi represents the i'th letter in the motif
wolffd@0: % and B is the background state.
wolffd@0: % Si chooses between running the motif or the background.
wolffd@0: % The Si and B states have self loops (not shown).
wolffd@0: %
wolffd@0: % The transition params are defined to respect the above topology.
wolffd@0: % The background is uniform; each motif state has a random obs. distribution.
wolffd@0: %
wolffd@0: % BNET = MK_MOTIF_HHMM(MOTIF_LENGTH, MOTIF_PATTERN)
wolffd@0: % In this case, we make the motif submodel deterministically
wolffd@0: % emit the motif pattern. 
wolffd@0: %
wolffd@0: % BNET = MK_MOTIF_HHMM(MOTIF_LENGTH, MOTIF_PATTERN, BACKGROUND_CHAR)
wolffd@0: % In this case, we make the background submodel
wolffd@0: % deterministically emit the specified character (to make the pattern
wolffd@0: % easier to see).
wolffd@0: 
wolffd@0: if nargin < 2, motif_pattern = []; end
wolffd@0: if nargin < 3, background_char = []; end
wolffd@0: 
wolffd@0: chars = ['a', 'c', 'g', 't'];
wolffd@0: Osize = length(chars);
wolffd@0: 
wolffd@0: motif_length = length(motif_pattern);
wolffd@0: Qsize = [2 motif_length];
wolffd@0: Qnodes = 1:2;
wolffd@0: D = 2;
wolffd@0: transprob = cell(1,D);
wolffd@0: termprob = cell(1,D);
wolffd@0: startprob = cell(1,D);
wolffd@0: 
wolffd@0: % startprob{d}(k,j), startprob{1}(1,j)
wolffd@0: % transprob{d}(i,k,j), transprob{1}(i,j)
wolffd@0: % termprob{d}(k,j)
wolffd@0: 
wolffd@0: 
wolffd@0: % LEVEL 1
wolffd@0: 
wolffd@0: startprob{1} = zeros(1, 2);
wolffd@0: startprob{1} = [1 0]; % always start in the background model
wolffd@0: 
wolffd@0: % When in the background state, we stay there with high prob
wolffd@0: % When in the motif state, we immediately return to the background state.
wolffd@0: transprob{1} = [0.8 0.2;
wolffd@0: 		1.0 0.0];
wolffd@0: 
wolffd@0: 
wolffd@0: % LEVEL 2
wolffd@0: startprob{2} = 'leftstart'; % both submodels start in substate 1
wolffd@0: transprob{2} = zeros(motif_length, 2, motif_length);
wolffd@0: termprob{2} = zeros(2, motif_length);
wolffd@0: 
wolffd@0: % In the background model, we only use state 1.
wolffd@0: transprob{2}(1,1,1) = 1; % self loop
wolffd@0: termprob{2}(1,1) = 0.2; % prob transition to end state
wolffd@0: 
wolffd@0: % Motif model
wolffd@0: transprob{2}(:,2,:) = mk_leftright_transmat(motif_length, 0); % no self loops
wolffd@0: termprob{2}(2,end) = 1.0; % last state immediately terminates
wolffd@0: 
wolffd@0: 
wolffd@0: % OBS LEVEl
wolffd@0: 
wolffd@0: obsprob = zeros([Qsize Osize]);
wolffd@0: if isempty(background_char)
wolffd@0:   % uniform background model
wolffd@0:   obsprob(1,1,:) = normalise(ones(Osize,1));
wolffd@0: else
wolffd@0:   % deterministic background model (easy to see!)
wolffd@0:   m = find(chars==background_char);
wolffd@0:   obsprob(1,1,m) = 1.0;
wolffd@0: end
wolffd@0: 
wolffd@0: if gen_motif
wolffd@0:   % initialise with true motif (cheating)
wolffd@0:   for i=1:motif_length
wolffd@0:     m = find(chars == motif_pattern(i));
wolffd@0:     obsprob(2,i,m) = 1.0;
wolffd@0:   end
wolffd@0: else
wolffd@0:   obsprob(2,:,:) = mk_stochastic(ones(motif_length, Osize));
wolffd@0: end
wolffd@0: 
wolffd@0: Oargs = {'CPT', obsprob};
wolffd@0: 
wolffd@0: [bnet, Qnodes, Fnodes, Onode] = mk_hhmm('Qsizes', Qsize, 'Osize', Osize, 'discrete_obs', 1, ...
wolffd@0: 	       'Oargs', Oargs, 'Ops', Qnodes(1:2), ...
wolffd@0: 	       'startprob', startprob, 'transprob', transprob, 'termprob', termprob);
wolffd@0: