## 2/06/2013

### HMM forward algorithm - c source code

This is HMM forward algorithm source code.

I wouldn't introduce HMM and forward algorithm on this page. Because there are so many good material in the web.

The source code need Txt file included A, B, pi and observation values.
After load the txt file, the code return likelihood value.

Refer tp the real excute screen.

txt file include...
number of state, number of observation, number of observation sequence
state matrix by 1 line
observation matrix by 1 line
observation sequence

This picture describes the result of forward algorithm.

C source code

```#include < stdio.h>

double HMM_forwardAlgorithm(double* &pi, double** &A, double ** &B, int* &obS, int sN, int oN, int obs_N);

void main()
{
///////////////////////////////////////////////
//Input HMM ramda
int stateN, observeN, ob_seq_N;
double *pi;
double **A;
double **B;
int *obS;
////////////////////////////////////////////////

FILE * fp;
fp = fopen("HMM_Ramda2.txt", "r");
fscanf(fp,"%d %d %d",&stateN, &observeN, &ob_seq_N);
/////////////////////////////////////////////////
//alloc buffer
pi = new double[stateN];
A = new double*[stateN];
B = new double*[stateN];
obS = new int[ob_seq_N];

for(int i=0; i< stateN; ++i)
{
A[i] = new double[stateN];
B[i] = new double[observeN];
}
/////////////////////////////////////////////////
//pi
for(int i=0; i< stateN; ++i)
fscanf(fp, "%lf", &(pi[i]) );
//A
for(int i=0; i< stateN; ++i)
for(int j=0; j< stateN; ++j)
fscanf(fp, "%lf", &(A[i][j]) );
//B
for(int i=0; i< stateN; ++i)
for(int j=0; j< observeN; ++j)
fscanf(fp, "%lf", &(B[i][j]) );

//observe sequence
for(int i=0; i< ob_seq_N; ++i)
fscanf(fp, "%d", &(obS[i]));

///////////////////////////////////////////////////
HMM_forwardAlgorithm(pi, A, B, obS, stateN, observeN, ob_seq_N);

////////////////////////////////////////////////
fclose(fp);

for(int i=0; i< stateN; ++i)
{
delete[] A[i];
delete[] B[i];
}
delete[] A;
delete[] B;
delete[] pi;
delete[] obS;

}

double HMM_forwardAlgorithm(double* &pi, double** &A, double ** &B, int* &obS, int sN, int oN, int obs_N)
{
int i,j,t;
//make forward matrix
double** fMtx = new double *[sN];
for(i=0; i< sN; ++i)
fMtx[i] = new double[obs_N];

//first step
for(i=0; i< sN; ++i)
{
fMtx[i][0] = pi[i] * B[i][ obS[0]-1 ];
}

//routine
double sum;
for(t=1; t< obs_N; ++t)
{
for(j=0; j< sN; ++j)
{
sum=0;
for(i=0; i< sN; ++i)
{
sum += (fMtx[i][t-1]*A[i][j]);
}
fMtx[j][t] = sum * B[j][ obS[t]-1];
}
}

//final
sum=0;
for(j=0; j< sN; ++j)
{
sum += (fMtx[j][obs_N-1]);
}

////////////////////////////////////////////////////////
//report
printf("state N=%d / Observe N=%d\n\n", sN, oN);

printf("Initial Matrix\n");
for(int i=0; i< sN; ++i)
printf("%lf ", pi[i] );
printf("\n\n");

printf("A matrix\n");
for(int i=0; i< sN; ++i)
{
for(int j=0; j< sN; ++j)
{
printf("%lf ", A[i][j]);
}
printf("\n");
}
printf("\n");
printf("B matrix\n");
for(int i=0; i< sN; ++i)
{
for(int j=0; j< oN; ++j)
{
printf("%lf ", B[i][j]);
}
printf("\n");
}
printf("\n");
printf("Observation sequence\n");
for(int i=0; i< obs_N; ++i)
printf("%d ", obS[i]);

printf("\n\n");

printf("** Trellis matrix **\n");
for(j=0; j< sN; ++j)
{
for(t=0; t< obs_N; ++t)
{
printf("%lf  ", fMtx[j][t]);

if(t!=obs_N-1)
printf("->  ");
}
printf("\n");
}

printf("\nLikelihood = %lf \n", sum);

for(i=0; i< sN; ++i)
delete[] fMtx[i];
delete[] fMtx;

return sum;
}

```