markovchain.c,v

Markov 过程

head	1.1;
access;
symbols;
locks
	olethros:1.1; strict;
comment	@ * @;


1.1
date	2003.02.09.18.37.34;	author olethros;	state Exp;
branches;
next	;


desc
@Markov Chain main
@


1.1
log
@Initial revision
@
text
@/* -*- Mode: C++; tab-width: 2; indent-tabs-mode: nil; c-basic-offset: 2 -*- */
/*VER: $Id: SampleBuffer.h,v 1.9 2002/12/20 00:28:01 olethros Exp olethros $*/

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <assert.h>
#include "debug.h"
#include "MarkovChain.h"


//==========================================================
// NewMarkovChain()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Creates a new markov chain.

  ARGUMENTS

  - n_states: The distinct observations possible. For example, if we
    wish to model transitions between letters, n_states might be equal
    to the number of letters in the alphabet, plus the number of
    punctuation symbols we use to incorporate.

  - mem_size: How large the state history should be. If 1, then only
    the previous observation is used to calculate transition
    probabilities. If n, then transition probabilities are considered
    between n-tuples of observation and the next observation.

  RETURN VALUE

  - A pointer to the newly created chain, NULL if failed.

  SEE ALSO

  - DeleteMarkovChain()
*/
MarkovChain* NewMarkovChain (int n_states, int mem_size)
{
  MarkovChain* chain = NULL;
  int n_transitions = 0;
  int i;

  if ((chain = AllocM (MarkovChain, 1))==NULL) {
    return NULL;
  }

  if ((chain->memory = AllocM (int, n_states))==NULL) {
    goto quick_exit;
  }

  chain->n_states = n_states;
  chain->mem_size = mem_size;
  chain->tot_states = 1;

  for (i=0; i<mem_size; i++) {
    chain->tot_states *= n_states;
  }

  chain->normalized = false;

  printf ("External States: %d\n", n_states);
  printf ("MemSize: %d State size : %d\n", chain->mem_size, chain->tot_states);
  printf ("Number of Possible Transitions: %d\n", n_states*chain->tot_states);
  n_transitions = n_states*chain->tot_states;
  if ((chain->transitions = AllocM (float, n_transitions))==NULL) {
    goto quick_exit;
  }
  /* Clear up the memory before use! */
  for (i=0; i<mem_size; i++) {
    chain->memory[i] = 0;
  }

  for (i=0; i<n_transitions; i++) {
    chain->transitions[i] = 0.0f;
  }
    
  return chain;

 quick_exit:
  DeleteMarkovChain (chain);
  return NULL;
}


//==========================================================
// DeleteMarkovChain()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Frees everything related to the markov chain.

  ARGUMENTS

  - chain: A pointer to the chain that you wish to free.

  RETURN VALUE

  - 0 if everything was OK, otherwise -1.

*/
int DeleteMarkovChain (MarkovChain* chain) {

  if (chain) {

    if (chain->transitions) {
      FreeM (chain->transitions);
    } else {
      fprintf (stderr, "Chain->transitions pointer null\n");
      return -1;
    }

    if (chain->memory) {
      FreeM (chain->memory);
    } else {
      fprintf (stderr, "Chain->memory pointer null\n");
      return -1;
    }
    
    FreeM (chain);

  } else {
    fprintf (stderr, "Chain pointer null\n");
    return -1;
  }

  return 0;
}



//==========================================================
// MarkovChainTrain()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Incrementally trains the markov chain with a new
    observations. When wishing to add the statistics of a particular
    sequence to the transition probabilities, first call
    MarkovChainReset() and then call MarkovChainTrain() iteratively,
    presenting elements of the sequence in order.

  ARGUMENTS

  - chain: A pointer to the chain that you wish to train.
  - state: The state ID.

  RETURN VALUE

  - Correct operation is indicated by 0.

  SEE ALSO

  - MarkovChainNormalize(), MarkovChainReset()

*/
int MarkovChainTrain (MarkovChain* chain, int state)
{
  int curr_state;
  
  assert((state>=0)&&(state<chain->n_states));

  curr_state = MarkovChainCalculateStateID (chain);
  MarkovChainPushState (chain, state);
  chain->transitions[state * chain->tot_states + curr_state] += 1.0f;

  return 0;
}


//==========================================================
// MarkovChainNormalize()
//----------------------------------------------------------
/*
  DESCRIPTION

  - After the statistic of all the sequences have been collected
    through the use of MarkovChainTrain(), call this function to
    convert the statistics to transition probabilities. This must be
    done before MarkovChainGenerate() is called.

  ARGUMENTS

  - 

  RETURN VALUE

  - 

  SEE ALSO

  - 
*/
int MarkovChainNormalize (MarkovChain* chain)
{
  // To speed mem-access a bit (also less typing :P)
  int n_states = chain->tot_states;
  float* trans = chain->transitions;
  int i;
  int spec=0;

  for (i=0; i<chain->tot_states; i++) {
    float tot =0.0f;
    int j;
    
    for (j=0; j<chain->n_states; j++) {
      if (trans[i+j*n_states]>0.0f) {
        spec++;
      }
      tot += trans[i+j*n_states];
    }
   
    if (tot>0.0f) {
      for (j=0; j<chain->n_states; j++) {
        trans[i+j*n_states] = trans[i+j*n_states]/tot;
      }
    } else {
      for (j=0; j<chain->n_states; j++) {
        trans[i+j*n_states] = 1.0f/((float) (chain->n_states));
      }
    }
  }

  chain->normalized = true;
  return spec;
}

//==========================================================
// MarkovChainReset()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Sets the state (and history) to 0. Call before training for
    distinct sequences and before generating a new
    sequence. (Otherwise training and generation will depend on past
    sequences).

  ARGUMENTS

  - chain: a pointer to the markov chain.

  RETURN VALUE

  - 0 for succesfull completion
*/
int MarkovChainReset (MarkovChain* chain) {
  int i;
  for (i=0; i<chain->mem_size; i++) {
    chain->memory[i] = 0;
  }
  chain->curr_addr = 0;
  return 0;
}

//==========================================================
// MarkovChainGenerate()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Generates a new observation based on the current state history,
  according to the transition tables. Note that you might generate as
  many new states as you wish. The state history must be explicitly
  updated by calling MarkovChainPushState() with the generated value
  (or one of the generated values, if you have called this multiple
  times, or if you are selecting a generated value from a number of
  different markov chains).

  ARGUMENTS

  - chain: A pointer to the chain

  RETURN VALUE

  - The ID of the next state, as generated by the MC. Returns -1 if
    nothing could be generated.

*/
int MarkovChainGenerate (MarkovChain* chain) {
  float* trans = chain->transitions;
  int n_states = chain->tot_states;
  int j;
  int curr;
  float sel;
  float tot = 0.0f;

  if (chain->normalized==false) {
    fprintf (stderr, "Warning: The chain was normalized (This error will not appear again)\n");
    chain->normalized = true;
  }

  curr = MarkovChainCalculateStateID (chain);
  sel = drand48();
  for (j=0; j<chain->n_states; j++) {
    //    printf ("(%f)",trans[curr+j*n_states]);
    tot += trans[curr+j*n_states];
    if ((sel<=tot)&&(trans[curr+j*n_states]>0.0f)) {

      return j;
    }
  }

  return -1;
}

//==========================================================
// MarkovChainCalculateStateID()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Calculates a unique state ID from the MC's state history.
  The calculation is $\sum_i s_{t-i} N^i$, where $N$ is the number of
  states, $i \in [0,M-1]$, where $M$ is the history size (the order of
  the markov chain) and $s_{t}$ is the state at time $t$.

  ARGUMENTS

  - chain: An argument

  RETURN VALUE

  - The id of the current state history.

  SEE ALSO

  - MarkovChainPushState(), MarkovChainReset()
*/
int MarkovChainCalculateStateID (MarkovChain* chain) {
  int id = 0;
  int i;
  int n = 1;

  //    printf ("ID[");
  for (i=1; i<=chain->mem_size; i++, n*=chain->n_states) {
    //    printf ("%d,",chain->memory[i-1]);
    id += chain->memory[i-1]*n; 
  }
  //        printf("]= %d\n",id);
  return id;
}


//==========================================================
// MarkovChainPushState()
//----------------------------------------------------------
/*
  DESCRIPTION

  - Adds a new state to the history FIFO, pushing out the oldest
    member of the history. Note that the implementation could be
    faster, but the history is pretty small. If you change the
    implementation you also need to change
    MarkovChainCalculateStateID() at least.

  ARGUMENTS

  - chain: A pointer to the chain.
  - state: A new state

  RETURN VALUE

  - It returns the ID of the popped state.
*/
int  MarkovChainPushState (MarkovChain* chain, int state) {
  int i;
  int popped_state;

  popped_state = chain->memory[chain->mem_size - 1];
  for (i=chain->mem_size - 1; i>0; i--) {
    chain->memory[i] = chain->memory[i-1];
  }
  chain->memory[0] = state;
  return popped_state;
}

//==========================================================
// MarkovChainShowTransitions)
//----------------------------------------------------------
/*
  DESCRIPTION

  - Can be useful for debugging.
*/

int MarkovChainShowTransitions (MarkovChain* chain) {
  float* trans = chain->transitions;
  int n_states = chain->tot_states;
  int j;
  int curr;
  float tot = 0.0f;

  printf ("\nState transition dump\n\n");
  curr = MarkovChainCalculateStateID (chain);
  printf ("Current state ID : %d\n", curr);
  printf ("Transition probabilities for next states:\n");

  for (j=0; j<chain->n_states; j++) {
    printf ("->%d = %f\n", j, trans[curr+j*n_states]);
    tot += trans[curr+j*n_states];
  }
  printf ("\nTOTAL: %f\n\n",tot);
  return 0;
}

@