Skip to content

AFNI/NIfTI Server

Sections
Personal tools
You are here: Home » AFNI » Documentation

Doxygen Source Code Documentation


Main Page   Alphabetical List   Data Structures   File List   Data Fields   Globals   Search  

Simplexx.c File Reference

Go to the source code of this file.


Functions

void allocate_arrays (float ***simplex, float **centroid, float **response, float **step_size, float **test1, float **test2)
void deallocate_arrays (float ***simplex, float **centroid, float **response, float **step_size, float **test1, float **test2)
void eval_vertices (float *response, int *worst, int *next, int *best)
void restart (float **simplex, float *response, float *step_size)
void calc_centroid (float **simplex, int worst, float *centroid)
void calc_reflection (float **simplex, float *centroid, int worst, float coef, float *vertex)
void replace (float **simplex, float *response, int index, float *vertex, float resp)
float calc_good_fit (float *response)
void simplex_initialize (float *parameters, float **simplex, float *response, float *step_size)
void simplex_optimization (float *parameters, float *sse)

Variables

int count = 0
int number_restarts = 0

Function Documentation

void allocate_arrays float ***    simplex,
float **    centroid,
float **    response,
float **    step_size,
float **    test1,
float **    test2
 

Definition at line 39 of file Simplexx.c.

References i, malloc, and test1().

Referenced by simplex_optimization().

00042 {
00043   int i;
00044 
00045   *centroid = (float *) malloc (sizeof(float) * DIMENSION);
00046   *response = (float *) malloc (sizeof(float) * (DIMENSION+1));
00047   *step_size = (float *) malloc (sizeof(float) * DIMENSION);
00048   *test1 = (float *) malloc (sizeof(float) * DIMENSION);
00049   *test2 = (float *) malloc (sizeof(float) * DIMENSION);
00050    
00051   *simplex = (float **) malloc (sizeof(float *) * (DIMENSION+1));
00052 
00053   for (i = 0;  i < DIMENSION+1;  i++)
00054     (*simplex)[i] = (float *) malloc (sizeof(float) * DIMENSION);
00055        
00056 }

void calc_centroid float **    simplex,
int    worst,
float *    centroid
 

Definition at line 157 of file Simplexx.c.

References i.

Referenced by simplex_optimization().

00158 {
00159   int i, j;
00160 
00161   for (i = 0;  i < DIMENSION;  i++)
00162     {
00163       centroid[i] = 0.0;
00164 
00165       /* add each vertex, except the worst */
00166       for (j = 0;  j < DIMENSION+1;  j++)
00167         if (j != worst)
00168           centroid[i] += simplex[j][i];
00169     }
00170 
00171   /* divide by the number of vertices */
00172   for (i = 0;  i < DIMENSION;  i++)
00173     centroid[i] /= DIMENSION;
00174 }

float calc_good_fit float *    response
 

Definition at line 214 of file Simplexx.c.

References avg, and i.

Referenced by simplex_optimization().

00215 {
00216   int i;
00217 
00218   float avg, sd, tmp;
00219 
00220   /* average the response values */
00221   avg = 0.0;
00222   for (i = 0;  i < DIMENSION+1;  i++)
00223     avg += response[i];
00224   avg /= DIMENSION+1;
00225 
00226   /* compute standard deviation of response */
00227   sd = 0.0;
00228   for (i = 0;  i < DIMENSION+1;  i++)
00229     {
00230       tmp = response[i] - avg;
00231       sd += tmp*tmp;
00232     }
00233   sd /= DIMENSION;
00234 
00235   return (sqrt(sd));
00236 }

void calc_reflection float **    simplex,
float *    centroid,
int    worst,
float    coef,
float *    vertex
 

Definition at line 182 of file Simplexx.c.

References i.

Referenced by simplex_optimization().

00184 {
00185   int i;
00186 
00187   for (i = 0;  i < DIMENSION;  i++)
00188     vertex[i] = centroid[i] + coef*(centroid[i] - simplex[worst][i]);
00189 }

void deallocate_arrays float ***    simplex,
float **    centroid,
float **    response,
float **    step_size,
float **    test1,
float **    test2
 

Definition at line 61 of file Simplexx.c.

References free, i, and test1().

Referenced by simplex_optimization().

00064 {
00065   int i;
00066 
00067   free (*centroid);    *centroid = NULL;
00068   free (*response);    *response = NULL;
00069   free (*step_size);   *step_size = NULL;
00070   free (*test1);       *test1 = NULL;
00071   free (*test2);       *test2 = NULL;
00072 
00073   for (i = 0;  i < DIMENSION+1;  i++)
00074     {
00075       free ((*simplex)[i]);
00076       (*simplex)[i] = NULL;
00077     }
00078 
00079   free (*simplex);     *simplex = NULL;
00080        
00081 }

void eval_vertices float *    response,
int *    worst,
int *    next,
int *    best
 

Definition at line 86 of file Simplexx.c.

References i.

Referenced by restart(), simplex_initialize(), and simplex_optimization().

00087 {
00088   int i;
00089 
00090   /* initialize values */
00091   *worst = 0;
00092   *best = 0;
00093 
00094   /* find the best and worst */
00095   for (i = 1;  i < DIMENSION+1;  i++)
00096     {
00097       if (response[i] > response[*worst])
00098         *worst = i;
00099       if (response[i] < response[*best])
00100         *best = i;
00101     }
00102 
00103   /* find the next worst index */
00104   if (*worst == 0)
00105     *next = 1;
00106   else
00107     *next = 0;
00108   
00109   for (i = 0;  i < DIMENSION+1;  i++)
00110     if ((i != *worst) && (response[i] > response[*next]))
00111       *next = i;
00112 }

void replace float **    simplex,
float *    response,
int    index,
float *    vertex,
float    resp
 

Definition at line 197 of file Simplexx.c.

References i.

Referenced by simplex_initialize(), and simplex_optimization().

00199 {
00200   int i;
00201 
00202   for (i = 0;  i < DIMENSION;  i++)
00203     simplex[index][i] = vertex[i];
00204 
00205   response[index] = resp;
00206 }

void restart float **    simplex,
float *    response,
float *    step_size
 

Definition at line 117 of file Simplexx.c.

References calc_error(), eval_vertices(), i, maxval, and rand_uniform().

Referenced by simplex_optimization().

00119 {
00120   const float STEP_FACTOR = 0.9;
00121   int i, j;
00122   int worst, next, best;
00123   float minval, maxval;
00124 
00125 
00126   /* find the current best vertex */
00127   eval_vertices (response, &worst, &next, &best); 
00128 
00129   /* set the first vertex to the current best */
00130   for (i = 0; i < DIMENSION;  i++)
00131     simplex[0][i] = simplex[best][i];
00132 
00133   /* decrease step size */
00134   for (i = 0;  i < DIMENSION;  i++)
00135     step_size[i] *= STEP_FACTOR;
00136 
00137   /* set up remaining vertices of simplex using new step size */
00138   for (i = 1;  i < DIMENSION+1;  i++)
00139     for (j = 0;  j < DIMENSION;  j++)
00140       {
00141         minval = simplex[0][j] - step_size[j];
00142         maxval = simplex[0][j] + step_size[j];
00143       simplex[i][j] = rand_uniform (minval, maxval);
00144       }
00145 
00146   /* initialize response for each vector */
00147   for (i = 0;  i < DIMENSION+1;  i++)
00148     response[i] = calc_error (simplex[i]);
00149 }

void simplex_initialize float *    parameters,
float **    simplex,
float *    response,
float *    step_size
 

Definition at line 244 of file Simplexx.c.

References calc_error(), eval_vertices(), i, maxval, rand_uniform(), and replace().

Referenced by simplex_optimization().

00246 {
00247   int i, j;
00248   int worst, next, best;
00249   float resp;
00250   float minval, maxval;
00251 
00252 
00253   for (j = 0;  j < DIMENSION;  j++)
00254     {
00255       simplex[0][j] = parameters[j];
00256       step_size[j] = 0.5 * parameters[j];
00257     }
00258 
00259   for (i = 1;  i < DIMENSION+1;  i++)
00260     for (j = 0;  j < DIMENSION;  j++)
00261       {
00262         minval = simplex[0][j] - step_size[j];
00263         maxval = simplex[0][j] + step_size[j];
00264         simplex[i][j] = rand_uniform (minval, maxval);
00265       }
00266 
00267   for (i = 0;  i < DIMENSION+1;  i++)
00268       response[i] = calc_error (simplex[i]);
00269 
00270   for (i = 1;  i < 500;  i++)
00271     {
00272       for (j = 0;  j < DIMENSION;  j++)
00273         {
00274           minval = simplex[0][j] - step_size[j];
00275           maxval = simplex[0][j] + step_size[j];
00276           parameters[j] = rand_uniform (minval, maxval);
00277         }
00278 
00279       resp = calc_error (parameters);
00280       eval_vertices (response, &worst, &next, &best);
00281       if (resp < response[worst])
00282         replace (simplex, response, worst, parameters, resp);
00283     }
00284 }

void simplex_optimization float *    parameters,
float *    sse
 

Definition at line 294 of file Simplexx.c.

References allocate_arrays(), calc_centroid(), calc_error(), calc_good_fit(), calc_reflection(), deallocate_arrays(), eval_vertices(), fit, i, number_restarts, replace(), restart(), simplex_initialize(), and test1().

Referenced by calc_full_model(), estpdf(), estpdf_float(), and estpdf_short().

00295 {
00296   const int MAX_ITERATIONS = 100;
00297   const int MAX_RESTARTS = 25;
00298   const float EXPANSION_COEF = 2.0;
00299   const float REFLECTION_COEF = 1.0;
00300   const float CONTRACTION_COEF = 0.5;
00301   const float TOLERANCE = 1.0e-10;
00302 
00303   float ** simplex   = NULL;
00304   float * centroid   = NULL;
00305   float * response   = NULL;
00306   float * step_size  = NULL;
00307   float * test1      = NULL;
00308   float * test2      = NULL;
00309   float resp1, resp2;
00310   int i, worst, best, next;
00311   int num_iter, num_restarts;
00312   int done;
00313   float fit;
00314 
00315 
00316   allocate_arrays (&simplex, &centroid, &response, &step_size, &test1, &test2);
00317   
00318   simplex_initialize (parameters, simplex, response, step_size);
00319 
00320   /* start loop to do simplex optimization */
00321   num_iter = 0;
00322   num_restarts = 0;
00323   done = 0;
00324   
00325   while (!done)
00326     {
00327       /* find the worst vertex and compute centroid of remaining simplex, 
00328          discarding the worst vertex */
00329       eval_vertices (response, &worst, &next, &best);
00330       calc_centroid (simplex, worst, centroid);
00331       
00332       /* reflect the worst point through the centroid */
00333       calc_reflection (simplex, centroid, worst, 
00334                        REFLECTION_COEF, test1);
00335       resp1 = calc_error (test1);
00336 
00337       /* test the reflection against the best vertex and expand it if the
00338          reflection is better.  if not, keep the reflection */
00339       if (resp1 < response[best])
00340         {
00341           /* try expanding */
00342           calc_reflection (simplex, centroid, worst, EXPANSION_COEF, test2);
00343           resp2 = calc_error (test2);
00344           if (resp2 <= resp1)    /* keep expansion */     
00345             replace (simplex, response, worst, test2, resp2);
00346           else                   /* keep reflection */
00347             replace (simplex, response, worst, test1, resp1);
00348         }
00349       else if (resp1 < response[next])
00350         {
00351           /* new response is between the best and next worst 
00352              so keep reflection */
00353           replace (simplex, response, worst, test1, resp1); 
00354         }
00355           else
00356         {
00357           /* try contraction */
00358           if (resp1 >= response[worst])
00359             calc_reflection (simplex, centroid, worst, 
00360                              -CONTRACTION_COEF, test2);
00361           else
00362             calc_reflection (simplex, centroid, worst, 
00363                              CONTRACTION_COEF, test2);
00364           resp2 =  calc_error (test2);
00365           
00366           /* test the contracted response against the worst response */
00367           if (resp2 > response[worst])
00368             {
00369               /* new contracted response is worse, so decrease step size
00370                  and restart */
00371               num_iter = 0;
00372               num_restarts += 1;
00373               restart (simplex, response, step_size);
00374             }
00375           else       /* keep contraction */
00376             replace (simplex, response, worst, test2, resp2);
00377         }
00378 
00379       /* test to determine when to stop.  
00380          first, check the number of iterations */
00381       num_iter += 1;    /* increment iteration counter */
00382       if (num_iter >= MAX_ITERATIONS)
00383         {
00384           /* restart with smaller steps */
00385           num_iter = 0;
00386           num_restarts += 1;
00387           restart (simplex, response, step_size);
00388         }
00389 
00390       /* limit the number of restarts */
00391       if (num_restarts == MAX_RESTARTS)  done = 1;
00392 
00393       /* compare relative standard deviation of vertex responses 
00394          against a defined tolerance limit */
00395       fit = calc_good_fit (response);
00396       if (fit <= TOLERANCE)  done = 1;
00397 
00398       /* if done, copy the best solution to the output array */
00399       if (done) 
00400         {
00401           eval_vertices (response, &worst, &next, &best);
00402           for (i = 0;  i < DIMENSION;  i++)
00403             parameters[i] = simplex[best][i];
00404           *sse = response[best];
00405         }
00406 
00407     }  /* while (!done) */
00408  
00409   number_restarts = num_restarts;
00410   deallocate_arrays (&simplex, &centroid, &response, &step_size,
00411                      &test1, &test2);
00412 
00413 }

Variable Documentation

int count = 0
 

Definition at line 33 of file Simplexx.c.

int number_restarts = 0
 

Definition at line 34 of file Simplexx.c.

Referenced by simplex_optimization().

 

Powered by Plone

This site conforms to the following standards: