Doxygen Source Code Documentation
fim+.c File Reference
#include "RegAna.c"
#include "ranks.c"
Go to the source code of this file.
Defines | |
#define | MAX_OUTPUT_TYPE 12 |
#define | FIM_FitCoef (0) |
#define | FIM_BestIndex (1) |
#define | FIM_PrcntChange (2) |
#define | FIM_Baseline (4) |
#define | FIM_Correlation (6) |
#define | FIM_PrcntFromAve (3) |
#define | FIM_Average (5) |
#define | FIM_PrcntFromTop (7) |
#define | FIM_Topline (8) |
#define | FIM_SigmaResid (9) |
#define | FIM_SpearmanCC (10) |
#define | FIM_QuadrantCC (11) |
#define | USE_LEGENDRE |
Functions | |
double | legendre (double x, int m) |
int | init_indep_var_matrix (int q, int p, int NFirst, int N, int polort, int num_ort_files, int num_ideal_files, MRI_IMAGE **ort_array, int **ort_list, MRI_IMAGE **ideal_array, int **ideal_list, float *x_bot, float *x_ave, float *x_top, int *good_list, matrix *x) |
int | init_regression_analysis (int q, int p, int N, int num_idealts, matrix xdata, matrix *x_base, matrix *xtxinvxt_base, matrix *x_ideal, matrix *xtxinvxt_ideal, float **rarray) |
float | calc_CC (int N, float *x, float *y) |
float | calc_SpearmanCC (int N, float *r, float *s) |
float | sign (float x) |
float | calc_QuadrantCC (int N, float *r, float *s) |
float | percent_change (float num, float den) |
void | regression_analysis (int N, int q, int num_idealts, matrix x_base, matrix xtxinvxt_base, matrix *x_ideal, matrix *xtxinvxt_ideal, vector y, float *x_bot, float *x_ave, float *x_top, float **rarray, int *output_type, float *FimParams) |
void | report_results (int *output_type, float *FimParams, char **label) |
Variables | |
char * | OUTPUT_TYPE_name [MAX_OUTPUT_TYPE] |
char | lbuf [4096] |
char | sbuf [256] |
Define Documentation
|
Definition at line 48 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 45 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 43 of file fim+.c. Referenced by check_for_valid_inputs(), display_help_menu(), and regression_analysis(). |
|
Definition at line 46 of file fim+.c. Referenced by calculate_results(), display_help_menu(), regression_analysis(), and write_bucket_data(). |
|
Definition at line 42 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 44 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 47 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 49 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 53 of file fim+.c. Referenced by display_help_menu(), regression_analysis(), and write_bucket_data(). |
|
Definition at line 51 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 52 of file fim+.c. Referenced by display_help_menu(), regression_analysis(), and write_bucket_data(). |
|
Definition at line 50 of file fim+.c. Referenced by display_help_menu(), and regression_analysis(). |
|
Definition at line 34 of file fim+.c. Referenced by allocate_memory(), calculate_results(), get_options(), initialize_options(), report_results(), save_voxel(), terminate_program(), and write_bucket_data(). |
|
|
Function Documentation
|
Definition at line 462 of file fim+.c. References i. Referenced by calc_QuadrantCC(), and calc_SpearmanCC().
00468 { 00469 const float EPSILON = 1.0e-10; /* protect against divide by zero */ 00470 float cc; 00471 int i; 00472 float csum, xsum, ysum, prod; 00473 00474 00475 csum = xsum = ysum = 0.0; 00476 00477 for (i = 0; i < N; i++) 00478 { 00479 csum += x[i] * y[i]; 00480 xsum += x[i] * x[i]; 00481 ysum += y[i] * y[i]; 00482 } 00483 00484 00485 prod = xsum * ysum; 00486 if (prod < EPSILON) 00487 cc = 0.0; 00488 else 00489 cc = csum / sqrt(prod); 00490 00491 00492 return (cc); 00493 00494 } |
|
Definition at line 558 of file fim+.c. References calc_CC(), free, i, malloc, r, and sign. Referenced by regression_analysis().
00564 { 00565 float cc; 00566 float rbar; 00567 float * dr = NULL; 00568 float * ds = NULL; 00569 int i; 00570 00571 00572 dr = (float *) malloc (sizeof(float) * N); 00573 ds = (float *) malloc (sizeof(float) * N); 00574 00575 00576 rbar = (N+1.0) / 2.0; 00577 for (i = 0; i < N; i++) 00578 { 00579 dr[i] = sign(r[i] - rbar); 00580 ds[i] = sign(s[i] - rbar); 00581 } 00582 00583 00584 cc = calc_CC(N, dr, ds); 00585 00586 free (dr); dr = NULL; 00587 free (ds); ds = NULL; 00588 00589 return (cc); 00590 } |
|
Definition at line 503 of file fim+.c. References calc_CC(), free, i, malloc, and r. Referenced by regression_analysis().
00509 { 00510 float cc; 00511 float rbar; 00512 float * dr = NULL; 00513 float * ds = NULL; 00514 int i; 00515 00516 00517 dr = (float *) malloc (sizeof(float) * N); 00518 ds = (float *) malloc (sizeof(float) * N); 00519 00520 00521 rbar = (N+1.0) / 2.0; 00522 for (i = 0; i < N; i++) 00523 { 00524 dr[i] = r[i] - rbar; 00525 ds[i] = s[i] - rbar; 00526 } 00527 00528 00529 cc = calc_CC(N, dr, ds); 00530 00531 free (dr); dr = NULL; 00532 free (ds); ds = NULL; 00533 00534 return (cc); 00535 } |
|
if here, m > 20 ==> use recurrence relation * Definition at line 185 of file fim+.c. References BIG_NUMBER, far, i, legendre(), matrix_create(), matrix_destroy(), matrix_equate(), matrix_extract_rows(), matrix_initialize(), MRI_FLOAT_PTR, p, and q. Referenced by calculate_results().
00205 { 00206 const int BIG_NUMBER = 33333; 00207 int i; /* time index */ 00208 int m; /* X matrix column index */ 00209 int n; /* X matrix row index */ 00210 int is; /* input ideal index */ 00211 float * far = NULL; 00212 int nx, ny, iq, nq; 00213 int Ngood; 00214 matrix xgood; 00215 #ifdef USE_LEGENDRE 00216 double nfac,nsub ; 00217 #endif 00218 00219 00220 /*----- Initialize X matrix -----*/ 00221 matrix_create (N, p, x); 00222 matrix_initialize (&xgood); 00223 00224 #ifdef USE_LEGENDRE 00225 nsub = 0.5*(N-1) ; 00226 nfac = 1.0/nsub ; 00227 #endif 00228 00229 /*----- Set up columns of X matrix corresponding to polynomial orts -----*/ 00230 for (m = 0; m <= polort; m++) 00231 for (n = 0; n < N; n++) 00232 { 00233 #ifndef USE_LEGENDRE /** the old polort way: t^m **/ 00234 i = n + NFirst; 00235 (*x).elts[n][m] = pow ((double)i, (double)m); 00236 00237 #else /** the new polort way: Legendre P_m(t) - 29 Mar 2005 **/ 00238 00239 (*x).elts[n][m] = legendre( nfac*(n-nsub) , m ) ; 00240 #endif 00241 } 00242 00243 00244 /*----- Set up columns of X matrix corresponding to ort time series -----*/ 00245 for (is = 0; is < num_ort_files; is++) 00246 { 00247 far = MRI_FLOAT_PTR (ort_array[is]); 00248 nx = ort_array[is]->nx; 00249 ny = ort_array[is]->ny; 00250 00251 if (ort_list[is] == NULL) 00252 for (iq = 0; iq < ny; iq++) 00253 { 00254 for (n = 0; n < N; n++) 00255 { 00256 i = n + NFirst; 00257 (*x).elts[n][m] = *(far + iq*nx + i); 00258 } 00259 m++; 00260 } 00261 else 00262 { 00263 nq = ort_list[is][0]; 00264 for (iq = 1; iq <= nq; iq++) 00265 { 00266 for (n = 0; n < N; n++) 00267 { 00268 i = n + NFirst; 00269 (*x).elts[n][m] = *(far + ort_list[is][iq]*nx + i); 00270 } 00271 m++; 00272 } 00273 } 00274 } 00275 00276 00277 /*----- Set up columns of X matrix corresponding to ideal time series -----*/ 00278 for (is = 0; is < num_ideal_files; is++) 00279 { 00280 far = MRI_FLOAT_PTR (ideal_array[is]); 00281 nx = ideal_array[is]->nx; 00282 ny = ideal_array[is]->ny; 00283 00284 if (ideal_list[is] == NULL) 00285 for (iq = 0; iq < ny; iq++) 00286 { 00287 for (n = 0; n < N; n++) 00288 { 00289 i = n + NFirst; 00290 (*x).elts[n][m] = *(far + iq*nx + i); 00291 } 00292 00293 m++; 00294 } 00295 else 00296 { 00297 nq = ideal_list[is][0]; 00298 for (iq = 1; iq <= nq; iq++) 00299 { 00300 for (n = 0; n < N; n++) 00301 { 00302 i = n + NFirst; 00303 (*x).elts[n][m] = *(far + ideal_list[is][iq]*nx + i); 00304 } 00305 00306 m++; 00307 } 00308 } 00309 } 00310 00311 00312 /*----- Remove row if ideal contains value over 33333 -----*/ 00313 Ngood = N; 00314 m = 0; 00315 for (n = 0; n < N; n++) 00316 { 00317 for (is = q; is < p; is++) 00318 { 00319 if ((*x).elts[n][is] >= BIG_NUMBER) break; 00320 } 00321 if (is < p) 00322 { 00323 Ngood--; 00324 } 00325 else 00326 { 00327 good_list[m] = n; 00328 m++; 00329 } 00330 } 00331 matrix_extract_rows ((*x), Ngood, good_list, &xgood); 00332 matrix_equate (xgood, x); 00333 00334 00335 /*----- Find min, max, and ave for each column of the X matrix -----*/ 00336 for (is = 0; is < p; is++) 00337 { 00338 x_bot[is] = x_top[is] = (*x).elts[0][is]; 00339 x_ave[is] = 0.0; 00340 for (n = 0; n < Ngood; n++) 00341 { 00342 if ((*x).elts[n][is] < x_bot[is]) x_bot[is] = (*x).elts[n][is]; 00343 if ((*x).elts[n][is] > x_top[is]) x_top[is] = (*x).elts[n][is]; 00344 x_ave[is] += (*x).elts[n][is] / Ngood; 00345 } 00346 } 00347 00348 00349 matrix_destroy (&xgood); 00350 00351 return (Ngood); 00352 00353 } |
|
the new polort way: Legendre P_m(t) - 29 Mar 2005 * Definition at line 362 of file fim+.c. References calc_coef(), calc_matrices(), calc_resids(), column_to_vector(), vector::elts, free, malloc, matrix_destroy(), matrix_initialize(), MTEST, p, q, rank_darray(), vector_destroy(), and vector_initialize(). Referenced by calculate_results().
00370 : (1/(X'X))X' for baseline model */ 00371 matrix * x_ideal, /* extracted X matrices for ideal models */ 00372 matrix * xtxinvxt_ideal, /* matrix: (1/(X'X))X' for ideal models */ 00373 float ** rarray /* ranked arrays of ideal time series */ 00374 ) 00375 00376 { 00377 int * plist = NULL; /* list of model parameters */ 00378 int ip, it; /* parameter indices */ 00379 int is, js; /* ideal indices */ 00380 int jm; /* lag index */ 00381 int ok; /* flag for successful matrix calculation */ 00382 matrix xtxinv_temp; /* intermediate results */ 00383 vector ideal; /* ideal vector */ 00384 vector coef_temp; /* intermediate results */ 00385 vector xres; /* vector of residuals */ 00386 float sse_base; /* baseline model error sum of squares */ 00387 00388 00389 /*----- Initialize matrix -----*/ 00390 matrix_initialize (&xtxinv_temp); 00391 vector_initialize (&ideal); 00392 vector_initialize (&coef_temp); 00393 vector_initialize (&xres); 00394 00395 00396 /*----- Allocate memory -----*/ 00397 plist = (int *) malloc (sizeof(int) * p); MTEST (plist); 00398 00399 00400 /*----- Initialize matrices for the baseline model -----*/ 00401 for (ip = 0; ip < q; ip++) 00402 plist[ip] = ip; 00403 ok = calc_matrices (xdata, q, plist, x_base, &xtxinv_temp, xtxinvxt_base); 00404 if (!ok) { matrix_destroy (&xtxinv_temp); return (0); }; 00405 00406 00407 /*----- Initialize matrices for ideal functions -----*/ 00408 for (is = 0; is < num_idealts; is++) 00409 { 00410 for (ip = 0; ip < q; ip++) 00411 { 00412 plist[ip] = ip; 00413 } 00414 00415 plist[q] = q + is; 00416 00417 ok = calc_matrices (xdata, q+1, plist, 00418 &(x_ideal[is]), &xtxinv_temp, &(xtxinvxt_ideal[is])); 00419 if (!ok) { matrix_destroy (&xtxinv_temp); return (0); }; 00420 } 00421 00422 00423 /*----- Set up the ranked array for each ideal -----*/ 00424 for (is = 0; is < num_idealts; is++) 00425 { 00426 /*----- Convert column of matrix to vector -----*/ 00427 column_to_vector (xdata, q+is, &ideal); 00428 00429 /*----- Calculate regression coefficients for baseline model -----*/ 00430 calc_coef (*xtxinvxt_base, ideal, &coef_temp); 00431 00432 /*----- Calculate the error sum of squares for the baseline model -----*/ 00433 sse_base = calc_resids (*x_base, coef_temp, ideal, &xres); 00434 00435 /*----- Form rank array from residual array -----*/ 00436 rarray[is] = rank_darray (N, xres.elts); 00437 00438 } 00439 00440 00441 /*----- Destroy matrix -----*/ 00442 matrix_destroy (&xtxinv_temp); 00443 vector_destroy (&ideal); 00444 vector_destroy (&coef_temp); 00445 vector_destroy (&xres); 00446 00447 00448 /*----- Deallocate memory -----*/ 00449 free (plist); plist = NULL; 00450 00451 00452 return (1); 00453 } |
|
Definition at line 60 of file fim+.c. Referenced by init_indep_var_matrix(), and legendre().
00061 { 00062 if( m < 0 ) return 1.0 ; /* bad input */ 00063 00064 switch( m ){ /*** P_m(x) for m=0..20 ***/ 00065 case 0: return 1.0 ; 00066 case 1: return x ; 00067 case 2: return (3.0*x*x-1.0)/2.0 ; 00068 case 3: return (5.0*x*x-3.0)*x/2.0 ; 00069 case 4: return ((35.0*x*x-30.0)*x*x+3.0)/8.0 ; 00070 case 5: return ((63.0*x*x-70.0)*x*x+15.0)*x/8.0 ; 00071 case 6: return (((231.0*x*x-315.0)*x*x+105.0)*x*x-5.0)/16.0 ; 00072 case 7: return (((429.0*x*x-693.0)*x*x+315.0)*x*x-35.0)*x/16.0 ; 00073 case 8: return ((((6435.0*x*x-12012.0)*x*x+6930.0)*x*x-1260.0)*x*x+35.0)/128.0; 00074 00075 /** 07 Feb 2005: this part generated by Maple, then hand massaged **/ 00076 00077 case 9: 00078 return (0.24609375e1 + (-0.3609375e2 + (0.140765625e3 + (-0.20109375e3 00079 + 0.949609375e2 * x * x) * x * x) * x * x) * x * x) * x; 00080 00081 case 10: 00082 return -0.24609375e0 + (0.1353515625e2 + (-0.1173046875e3 + 00083 (0.3519140625e3 + (-0.42732421875e3 + 0.18042578125e3 * x * x) 00084 * x * x) * x * x) * x * x) * x * x; 00085 00086 case 11: 00087 return (-0.270703125e1 + (0.5865234375e2 + (-0.3519140625e3 + 00088 (0.8546484375e3 + (-0.90212890625e3 + 0.34444921875e3 * x * x) 00089 * x * x) * x * x) * x * x) * x * x) * x; 00090 00091 case 12: 00092 return 0.2255859375e0 + (-0.17595703125e2 + (0.2199462890625e3 + 00093 (-0.99708984375e3 + (0.20297900390625e4 + (-0.1894470703125e4 00094 + 0.6601943359375e3 * x * x) * x * x) * x * x) * x * x) * x * x) 00095 * x * x; 00096 00097 case 13: 00098 return (0.29326171875e1 + (-0.87978515625e2 + (0.7478173828125e3 + 00099 (-0.270638671875e4 + (0.47361767578125e4 + (-0.3961166015625e4 00100 + 0.12696044921875e4 * x * x) * x * x) * x * x) * x * x) * x * x) 00101 * x * x) * x; 00102 00103 case 14: 00104 return -0.20947265625e0 + (0.2199462890625e2 + (-0.37390869140625e3 + 00105 (0.236808837890625e4 + (-0.710426513671875e4 + 00106 (0.1089320654296875e5 + (-0.825242919921875e4 + 00107 0.244852294921875e4 * x * x) * x * x) * x * x) * x * x) * x * x) 00108 * x * x) * x * x; 00109 00110 case 15: 00111 return (-0.314208984375e1 + (0.12463623046875e3 + (-0.142085302734375e4 00112 + (0.710426513671875e4 + (-0.1815534423828125e5 + 00113 (0.2475728759765625e5 + (-0.1713966064453125e5 + 00114 0.473381103515625e4 * x * x) * x * x) * x * x) * x * x) 00115 * x * x) * x * x) * x * x) * x; 00116 00117 case 16: 00118 return 0.196380615234375e0 + (-0.26707763671875e2 + (0.5920220947265625e3 00119 + (-0.4972985595703125e4 + (0.2042476226806641e5 + 00120 (-0.4538836059570312e5 + (0.5570389709472656e5 + 00121 (-0.3550358276367188e5 + 0.9171758880615234e4 * x * x) * x * x) 00122 * x * x) * x * x) * x * x) * x * x) * x * x) * x * x; 00123 00124 case 17: 00125 return (0.3338470458984375e1 + (-0.169149169921875e3 + 00126 (0.2486492797851562e4 + (-0.1633980981445312e5 + 00127 (0.5673545074462891e5 + (-0.1114077941894531e6 + 00128 (0.1242625396728516e6 + (-0.7337407104492188e5 + 00129 0.1780400253295898e5 * x * x) * x * x) * x * x) * x * x) 00130 * x * x) * x * x) * x * x) * x * x) * x; 00131 00132 case 18: 00133 return -0.1854705810546875e0 + (0.3171546936035156e2 + 00134 (-0.8880331420898438e3 + (0.9531555725097656e4 + 00135 (-0.5106190567016602e5 + (0.153185717010498e6 + 00136 (-0.2692355026245117e6 + (0.275152766418457e6 + 00137 (-0.1513340215301514e6 + 0.3461889381408691e5 * x * x) * x * x) 00138 * x * x) * x * x) * x * x) * x * x) * x * x) * x * x) * x * x; 00139 00140 case 19: 00141 return (-0.3523941040039062e1 + (0.2220082855224609e3 + 00142 (-0.4084952453613281e4 + (0.3404127044677734e5 + 00143 (-0.153185717010498e6 + (0.4038532539367676e6 + 00144 (-0.6420231216430664e6 + (0.6053360861206055e6 + 00145 (-0.3115700443267822e6 + 0.6741574058532715e5 * x * x) * x * x) 00146 * x * x) * x * x) * x * x) * x * x) * x * x) * x * x) * x * x) * x; 00147 00148 case 20: 00149 return 0.1761970520019531e0 + (-0.3700138092041016e2 + 00150 (0.127654764175415e4 + (-0.1702063522338867e5 + 00151 (0.1148892877578735e6 + (-0.4442385793304443e6 + 00152 (0.1043287572669983e7 + (-0.1513340215301514e7 + 00153 (0.1324172688388824e7 + (-0.6404495355606079e6 + 00154 0.1314606941413879e6 * x * x) * x * x) * x * x) * x * x) * x * x) 00155 * x * x) * x * x) * x * x) * x * x) * x * x; 00156 } 00157 00158 #if 0 00159 /* order out of range: return Chebyshev instead (it's easy) */ 00160 00161 if( x >= 1.0 ) x = 0.0 ; 00162 else if ( x <= -1.0 ) x = 3.14159265358979323846 ; 00163 else x = acos(x) ; 00164 return cos(m*x) ; 00165 #else 00166 /** if here, m > 20 ==> use recurrence relation **/ 00167 00168 { double pk, pkm1, pkm2 ; int k ; 00169 pkm2 = legendre( x , 19 ) ; 00170 pkm1 = legendre( x , 20 ) ; 00171 for( k=21 ; k <= m ; k++ , pkm2=pkm1 , pkm1=pk ) 00172 pk = ((2.0*k-1.0)*x*pkm1 - (k-1.0)*pkm2)/k ; 00173 return pk ; 00174 } 00175 #endif 00176 } |
|
Definition at line 599 of file fim+.c. References sign. Referenced by regression_analysis().
00600 { 00601 const float EPSILON = 1.0e-10; /* guard against divide by zero */ 00602 const float MAX_PERCENT = 1000.0; /* limit maximum percent change */ 00603 float PrcntChange; 00604 00605 00606 if (fabs(den) < EPSILON) 00607 PrcntChange = sign(num) * MAX_PERCENT; 00608 else 00609 PrcntChange = 100.0 * num / den; 00610 00611 if (PrcntChange > MAX_PERCENT) PrcntChange = MAX_PERCENT; 00612 if (PrcntChange < -MAX_PERCENT) PrcntChange = -MAX_PERCENT; 00613 00614 return (PrcntChange); 00615 } |
|
Definition at line 624 of file fim+.c. References calc_coef(), calc_QuadrantCC(), calc_resids(), calc_rsqr(), calc_SpearmanCC(), calc_sse(), vector::elts, FIM_Average, FIM_Baseline, FIM_BestIndex, FIM_Correlation, FIM_FitCoef, FIM_PrcntChange, FIM_PrcntFromAve, FIM_PrcntFromTop, FIM_QuadrantCC, FIM_SigmaResid, FIM_SpearmanCC, FIM_Topline, free, percent_change(), q, rank_darray(), vector_destroy(), vector_equate(), and vector_initialize(). Referenced by calculate_results().
00629 : (1/(X'X))X' for baseline model */ 00630 matrix * x_ideal, /* extracted X matrices for ideal models */ 00631 matrix * xtxinvxt_ideal, /* matrix: (1/(X'X))X' for ideal models */ 00632 vector y, /* vector of measured data */ 00633 float * x_bot, /* minimum of stimulus time series */ 00634 float * x_ave, /* average of stimulus time series */ 00635 float * x_top, /* maximum of stimulus time series */ 00636 float ** rarray, /* ranked arrays of ideal time series */ 00637 int * output_type, /* list of operator requested outputs */ 00638 float * FimParams /* output fim parameters */ 00639 ) 00640 00641 { 00642 const float EPSILON = 1.0e-05; /* protection against divide by zero */ 00643 int is; /* input ideal index */ 00644 float sse_base; /* error sum of squares, baseline model */ 00645 float sse_ideal; /* error sum of squares, ideal model */ 00646 vector coef_temp; /* intermediate results */ 00647 vector coef_best; /* best results */ 00648 vector yres; /* vector of residuals */ 00649 float rtemp, rbest; /* best correlation coefficient */ 00650 float stemp, sbest; /* best Spearman correlation coefficient */ 00651 float qtemp, qbest; /* best quadrant correlation coefficient */ 00652 float mse; /* mean square error (sample variance) */ 00653 00654 float * sarray = NULL; /* ranked array of measurement data */ 00655 00656 /* fim output parameters */ 00657 float FitCoef = 0.0; 00658 int BestIndex = 0; 00659 float PrcntChange = 0.0; 00660 float Baseline = 0.0; 00661 float Correlation = 0.0; 00662 float PrcntFromAve = 0.0; 00663 float Average = 0.0; 00664 float PrcntFromTop = 0.0; 00665 float Topline = 0.0; 00666 float SigmaResid = 0.0; 00667 float SpearmanCC = 0.0; 00668 float QuadrantCC = 0.0; 00669 00670 00671 /*----- Initialization -----*/ 00672 vector_initialize (&coef_temp); 00673 vector_initialize (&coef_best); 00674 vector_initialize (&yres); 00675 00676 00677 /*----- Calculate regression coefficients for baseline model -----*/ 00678 calc_coef (xtxinvxt_base, y, &coef_temp); 00679 00680 00681 /*----- Calculate the error sum of squares for the baseline model -----*/ 00682 sse_base = calc_resids (x_base, coef_temp, y, &yres); 00683 00684 00685 /*----- Form rank array from y array -----*/ 00686 if (output_type[FIM_SpearmanCC] || output_type[FIM_QuadrantCC]) 00687 sarray = rank_darray (N, yres.elts); 00688 00689 00690 /*----- Determine the best ideal reference for this voxel -----*/ 00691 rbest = 0.0; sbest = 0.0; qbest = 0.0; BestIndex = -1; 00692 for (is = 0; is < num_idealts; is++) 00693 { 00694 00695 /*----- Calculate regression coefficients for ideal model -----*/ 00696 calc_coef (xtxinvxt_ideal[is], y, &coef_temp); 00697 00698 00699 /*----- Calculate the error sum of squares for the ideal model -----*/ 00700 sse_ideal = calc_sse (x_ideal[is], coef_temp, y); 00701 00702 00703 /*----- Calculate partial R^2 for this ideal -----*/ 00704 rtemp = calc_rsqr (sse_ideal, sse_base); 00705 00706 00707 if (rtemp >= rbest) 00708 { 00709 rbest = rtemp; 00710 BestIndex = is; 00711 vector_equate (coef_temp, &coef_best); 00712 if (num_idealts == 1) 00713 mse = sse_ideal / (N-q-1); 00714 else 00715 mse = sse_ideal / (N-q-2); 00716 } 00717 00718 00719 /*----- Calculate the Spearman rank correlation coefficient -----*/ 00720 if (output_type[FIM_SpearmanCC]) 00721 { 00722 stemp = calc_SpearmanCC (N, rarray[is], sarray); 00723 if (fabs(stemp) > fabs(sbest)) sbest = stemp; 00724 } 00725 00726 00727 /*----- Calculate the Quadrant correlation coefficient -----*/ 00728 if (output_type[FIM_QuadrantCC]) 00729 { 00730 qtemp = calc_QuadrantCC (N, rarray[is], sarray); 00731 if (fabs(qtemp) > fabs(qbest)) qbest = qtemp; 00732 } 00733 } 00734 00735 if ((0 <= BestIndex) && (BestIndex < num_idealts)) 00736 { 00737 float Top, Ave, Base, Center; 00738 int ip; 00739 00740 Top = x_top[q+BestIndex]; 00741 Ave = x_ave[q+BestIndex]; 00742 Base = x_bot[q+BestIndex]; 00743 00744 00745 FitCoef = coef_best.elts[q]; 00746 Correlation = sqrt(rbest); 00747 if (FitCoef < 0.0) Correlation = -Correlation; 00748 00749 Center = 0.0; 00750 for (ip = 0; ip < q; ip++) 00751 Center += coef_best.elts[ip] * x_ave[ip]; 00752 00753 Baseline = Center + FitCoef*Base; 00754 Average = Center + FitCoef*Ave; 00755 Topline = Center + FitCoef*Top; 00756 00757 00758 00759 PrcntChange = percent_change (FitCoef * (Top-Base), Baseline); 00760 PrcntFromAve = percent_change (FitCoef * (Top-Base), Average); 00761 PrcntFromTop = percent_change (FitCoef * (Top-Base), Topline); 00762 00763 SigmaResid = sqrt(mse); 00764 00765 SpearmanCC = sbest; 00766 00767 QuadrantCC = qbest; 00768 } 00769 00770 00771 /*----- Save output parameters -----*/ 00772 FimParams[FIM_FitCoef] = FitCoef; 00773 FimParams[FIM_BestIndex] = BestIndex + 1.0; 00774 FimParams[FIM_PrcntChange] = PrcntChange; 00775 FimParams[FIM_Baseline] = Baseline; 00776 FimParams[FIM_Correlation] = Correlation; 00777 FimParams[FIM_PrcntFromAve] = PrcntFromAve; 00778 FimParams[FIM_Average] = Average; 00779 FimParams[FIM_PrcntFromTop] = PrcntFromTop; 00780 FimParams[FIM_Topline] = Topline; 00781 FimParams[FIM_SigmaResid] = SigmaResid; 00782 FimParams[FIM_SpearmanCC] = SpearmanCC; 00783 FimParams[FIM_QuadrantCC] = QuadrantCC; 00784 00785 00786 /*----- Dispose of vectors -----*/ 00787 vector_destroy (&coef_temp); 00788 vector_destroy (&coef_best); 00789 vector_destroy (&yres); 00790 00791 00792 /*----- Deallocate memory -----*/ 00793 if (sarray != NULL) 00794 { free (sarray); sarray = NULL; } 00795 00796 } |
|
Definition at line 806 of file fim+.c. References lbuf, MAX_OUTPUT_TYPE, OUTPUT_TYPE_name, and sbuf.
00812 { 00813 int ip; /* parameter index */ 00814 00815 00816 if( label != NULL ){ /* assemble this 1 line at a time from sbuf */ 00817 00818 lbuf[0] = '\0' ; /* make this a 0 length string to start */ 00819 00820 /** for each reference, make a string into sbuf **/ 00821 00822 for (ip = 0; ip < MAX_OUTPUT_TYPE; ip++) 00823 if (output_type[ip]) 00824 { 00825 sprintf (sbuf, "%12s = %10.4f \n", 00826 OUTPUT_TYPE_name[ip], FimParams[ip]); 00827 strcat (lbuf, sbuf); 00828 } 00829 00830 00831 *label = lbuf ; /* send address of lbuf back in what label points to */ 00832 } 00833 00834 } |
|
Definition at line 543 of file fim+.c.
00545 { 00546 if (x > 0.0) return (1.0); 00547 if (x < 0.0) return (-1.0); 00548 return (0.0); 00549 } |
Variable Documentation
|
Definition at line 801 of file fim+.c. Referenced by report_results(). |
|
Initial value: { "Fit Coef", "Best Index", "% Change", "% From Ave", "Baseline", "Average", "Correlation", "% From Top", "Topline", "Sigma Resid", "Spearman CC", "Quadrant CC" } Definition at line 36 of file fim+.c. Referenced by report_results(). |
|
Definition at line 802 of file fim+.c. Referenced by report_results(). |