fcm.cc

KNN algorithm, and related

/*--------------------------------------------------------------------------
@COPYRIGHT  :
              Copyright 1997, Vasken Kollokian, 
              McConnell Brain Imaging Centre,
              Montreal Neurological Institute, McGill University.
              Permission to use, copy, modify, and distribute this
              software and its documentation for any purpose and without
              fee is hereby granted, provided that the above copyright
              notice appear in all copies.  The author and McGill University
              make no representations about the suitability of this
              software for any purpose.  It is provided "as is" without
              express or implied warranty.
---------------------------------------------------------------------------- 
$RCSfile: fcm.cc,v $
$Revision: 1.1.1.1 $
$Author: jason $
$Date: 2002/03/20 22:16:34 $
$State: Exp $
--------------------------------------------------------------------------*/
/* ----------------------------- MNI Header -----------------------------------
@NAME       : fcm.c
@INPUT      : 
@OUTPUT     : 
@RETURNS    : 
@DESCRIPTION: fuzzy c means unsupervised classifier
@METHOD     : 
@GLOBALS    : 
@CALLS      : 
@CREATED    : Feb. 12, 1996 (Vasco KOLLOKIAN)
@MODIFIED   : $Log: fcm.cc,v $
@MODIFIED   : Revision 1.1.1.1  2002/03/20 22:16:34  jason
@MODIFIED   : first autoconfiscated version that compiles under linux gcc 3
@MODIFIED   :
@MODIFIED   : Revision 1.1.1.1  1997/02/11 00:06:48  alex
@MODIFIED   : Sources for classify, copied from Vasken Kollokian
@MODIFIED   :
 * Revision 1.3  1996/08/29  03:54:18  vasco
 * I forgot.
 *
 * Revision 1.2  1996/03/14  06:31:30  vasco
 * added static on some multiply defined variables.
 * removed some unused variables.
 *
 * Revision 1.1  1996/03/06  00:28:39  vasco
 * Initial revision
 *
---------------------------------------------------------------------------- */
extern "C" {
#include <volume_io.h>
#include <limits.h>
}
#include "../unsuper_globals.h" 
#include <unistd.h>
#include <sys/types.h>

/* locally defined prototypes */
void fcm_count_classes_and_set_names(void);

/* locally defined global variables */
static Real   **mean_feature_matrix;  /* matrix to reflect mean features */
static int    *mean_feature_class_vector;/* vector to reflect mean feature classes */
static Real   *euclidian_vector;         /* eucledian vector of each sample */
Real   *fuzzy_fcm_vector;               /* fuzzy vector of each sample */
static int    v1, v2, v3;               /* voxel pointers in x, y, and z direction */
static Real   m;                        /* proxemity neighbourhood index */

/* ----------------------------- MNI Header -----------------------------------
@NAME       : fcm_init_training
@INPUT      : 
@OUTPUT     : 
@RETURNS    : 
@DESCRIPTION: 
@METHOD     : 
@GLOBALS    : 
@CALLS      : 
@CREATED    : Feb 12, 1996 (Vasco KOLLOKIAN)
@MODIFIED   : 
---------------------------------------------------------------------------- */
void fcm_init_training(char *param_filename)
{


  FILE *fcm_par_file;

  /* if load_training is not used, then allocate space for struct */
  if ( !load_train_filename ) {

    /* count num_classes, and set their names */
    fcm_count_classes_and_set_names();

    /* reserve area for the euclidianvector */
    ALLOC( euclidian_vector, num_classes );

    /* reserve area for the fuzzy_fcm_vector */
    ALLOC( fuzzy_fcm_vector, num_classes );

    /* reserve area for the mean feature matrix */
    ALLOC2D(mean_feature_matrix, num_classes, num_features);

    /* reserve area for the mean_feature_class_vector */
    ALLOC(mean_feature_class_vector, num_classes);

  }

  /* check to see if the filename is there */
  if ( param_filename && !file_exists(param_filename)  ) {

    (void) fprintf(stderr,"File `%s' doesn't exist !\n ", param_filename);
    exit(EXIT_FAILURE);

  }

  if ( !param_filename ) {

    /* take the default here, in case a param file is not given */
    m = 2.0;

  }
  else {
    
    if (verbose) 
      fprintf(stdout, "Loading the parameter file %s\n", param_filename);      
    
    /* open the parameter file, and read the values  */
    fcm_par_file = fopen(param_filename, "r");
    
    if ( fcm_par_file == NULL) {

      fprintf(stderr, "Cannot open %s\n", param_filename);
      exit(EXIT_FAILURE);
    }
    
   /* scan for the neighbourhood proxemity number */
    fscanf( fcm_par_file, "m=%f\n", &m);
     
    fclose(fcm_par_file);

    if ( m <= 1.0 ) {
 
     fprintf(stderr, "Invalid m value : %f\n", m);
      exit(EXIT_FAILURE);
    }
 
  } /* else */

  if ( debug > 2) {

    fprintf(stdout, "m=%f\n", m);      
  }

}


/* ----------------------------- MNI Header -----------------------------------
@NAME       : fcm_train_samples
@INPUT      : 
@OUTPUT     : 
@RETURNS    : ?
@DESCRIPTION: takes a feature matrix and trains a classifier on it.
@METHOD     : 
@GLOBALS    : 
@CALLS      : 
@CREATED    : Feb 12, 1996 (Vasco KOLLOKIAN)
@MODIFIED   : 
---------------------------------------------------------------------------- */
void fcm_train_samples(void)
{
  
  int      i, j;              /* counters - samples, features, classes */  

  Real     u_ik, u_ikm, sig_u_ikm;


  /* initialize mean feature matrix, feature class vector & num of class samples */
  for_less( i, 0, num_classes) {

    for_less( j, 0, num_features) 
      mean_feature_matrix[i][j] = 0.0;

  }
   
 /* train on all input training volumes, one at a time, and one per class */
  for_less ( i, 0, num_classes ) {

    /* reset the u variable and the sig total */
    u_ik = 0.0;
    sig_u_ikm = 0.0;

    if (verbose)
      (void) fprintf(stdout, "Training on volume %s\n", trainvol_filename[i]);

    /* repeat for the total number of samples - the entire training volume */
    for_less( v1, 0, first_volume_sizes[0] ) {

      if ( verbose ) 
	write(2, "*", 1);

      for_less( v2, 0, first_volume_sizes[1] ) {

	if ( debug > 13 ) 
	  write(2, "#", 1);
	
	for_less( v3, 0, first_volume_sizes[2] ) {
	
	  /* get u_ik for class i, from the ith training volume */
	  u_ik =  get_volume_real_value(train_volume[i], 
					v1, 
					v2, 
					v3,
					0,0);
	
	  /* take it to the mth power  */
	  u_ikm = pow(u_ik, m);

	  /* calculate and store the total - Sig[ (uik)^m ] */
	  sig_u_ikm += u_ikm;

	  /* Sig [ (uik)^m * Xk ] for each and every feature */
	  for_less(j, 0, num_features)  
	    mean_feature_matrix[i][j] += u_ikm * get_volume_real_value(in_volume[j],
								       v1, 
								       v2, 
								       v3,
								       0,0);
	  
	    
	} /* for_less v3 */
      } /* for_less v2 */
    } /* for_less v1 */

    if (verbose)
      (void) fprintf(stdout, "\nGenerating class centroid for class %s\n",
		     class_name[i]);

    /* calculate for each Xk, Vi = (Sig[ (uik)^m *Xk ] ) / ( Sig[ (uik)^m ] ) */
    for_less(j, 0, num_features)  
      mean_feature_matrix[i][j] /= sig_u_ikm ;
    
    /* get rid of the training volume, once finished */
    delete_volume(train_volume[i]);

  } /* for_less ( i, 0, num_classes ) */

  if (debug > 2 ) {

    fprintf( stdout, "Printing class centroids  ...\n");

    for_less( i, 0, num_classes) {
      for_less( j, 0, num_features) 
	fprintf( stdout, "%f ", mean_feature_matrix[i][j]);
      fprintf( stdout, "%s\n", class_name[i]);      
    }

    fprintf( stdout, "-----------\n");

  }
  
} /* fcm_train_samples(void) */


/* ----------------------------- MNI Header -----------------------------------
@NAME       : fcm_load_training
@INPUT      : 
@OUTPUT     : 
@RETURNS    : 
@DESCRIPTION: 
@METHOD     : 
@GLOBALS    : 
@CALLS      : 
@CREATED    : Feb 12, 1996 (Vasco KOLLOKIAN)
@MODIFIED   : 
---------------------------------------------------------------------------- */
void fcm_load_training(char *load_train_filename)
{

  int i, j, k;