'Trains the HMM using the Baum-Welch algorithm to maximise the
probability of the data sequence. This is a variant of the EM
algorithm, and is unsupervised in that it doesn\'t need the state
sequences for the symbols. The code is based on \'A Tutorial on Hidden
Markov Models and Selected Applications in Speech Recognition\',
Lawrence Rabiner, IEEE, 1989.
:return: the trained model
:rtype: HiddenMarkovModelTagger
:param unlabeled_sequences: the training data, a set of
sequences of observations
:type unlabeled_sequences: list
kwargs may include following parameters:
:param model: a HiddenMarkovModelTagger instance used to begin
the Baum-Welch algorithm
:param max_iterations: the maximum number of EM iterations
:param convergence_logprob: the maximum change in log probability to
allow convergence'

def train_unsupervised(self, unlabeled_sequences, update_outputs=True, **kwargs):

    model = kwargs.get(u'model')
    if (not model):
        priors = RandomProbDist(self._states)
        transitions = DictionaryConditionalProbDist(dict(((state, RandomProbDist(self._states)) for state in self._states)))
        outputs = DictionaryConditionalProbDist(dict(((state, RandomProbDist(self._symbols)) for state in self._states)))
        model = HiddenMarkovModelTagger(self._symbols, self._states, transitions, outputs, priors)
    self._states = model._states
    self._symbols = model._symbols
    N = len(self._states)
    M = len(self._symbols)
    symbol_numbers = dict(((sym, i) for (i, sym) in enumerate(self._symbols)))
    model._transitions = DictionaryConditionalProbDist(dict(((s, MutableProbDist(model._transitions[s], self._states)) for s in self._states)))
    if update_outputs:
        model._outputs = DictionaryConditionalProbDist(dict(((s, MutableProbDist(model._outputs[s], self._symbols)) for s in self._states)))
    model.reset_cache()
    converged = False
    last_logprob = None
    iteration = 0
    max_iterations = kwargs.get(u'max_iterations', 1000)
    epsilon = kwargs.get(u'convergence_logprob', 1e-06)
    while ((not converged) and (iteration < max_iterations)):
        A_numer = _ninf_array((N, N))
        B_numer = _ninf_array((N, M))
        A_denom = _ninf_array(N)
        B_denom = _ninf_array(N)
        logprob = 0
        for sequence in unlabeled_sequences:
            sequence = list(sequence)
            if (not sequence):
                continue
            (lpk, seq_A_numer, seq_A_denom, seq_B_numer, seq_B_denom) = self._baum_welch_step(sequence, model, symbol_numbers)
            for i in range(N):
                A_numer[i] = np.logaddexp2(A_numer[i], (seq_A_numer[i] - lpk))
                B_numer[i] = np.logaddexp2(B_numer[i], (seq_B_numer[i] - lpk))
            A_denom = np.logaddexp2(A_denom, (seq_A_denom - lpk))
            B_denom = np.logaddexp2(B_denom, (seq_B_denom - lpk))
            logprob += lpk
        for i in range(N):
            logprob_Ai = (A_numer[i] - A_denom[i])
            logprob_Bi = (B_numer[i] - B_denom[i])
            logprob_Ai -= logsumexp2(logprob_Ai)
            logprob_Bi -= logsumexp2(logprob_Bi)
            si = self._states[i]
            for j in range(N):
                sj = self._states[j]
                model._transitions[si].update(sj, logprob_Ai[j])
            if update_outputs:
                for k in range(M):
                    ok = self._symbols[k]
                    model._outputs[si].update(ok, logprob_Bi[k])
        if ((iteration > 0) and (abs((logprob - last_logprob)) < epsilon)):
            converged = True
        print(u'iteration', iteration, u'logprob', logprob)
        iteration += 1
        last_logprob = logprob
    return model
