1 /***************************************************************************** 2 Major portions of this software are copyrighted by the Medical College 3 of Wisconsin, 1994-2000, and are released under the Gnu General Public 4 License, Version 2. See the file README.Copyright for details. 5 ******************************************************************************/ 6 7 /* 8 This is the header file for NLfit.c. 9 10 File: NLfit.h 11 Author: B. Douglas Ward 12 Date: 3 June 1997 13 14 Mod: Added options for percent signal change above baseline, and 15 calculation of area under the signal above baseline. 16 Date: 26 November 1997 17 18 Mod: Added novar flag to eliminate unnecessary calculations. 19 Date: 13 July 1999 20 21 Mod: Changes for output of R^2 (coefficient of multiple determination), 22 and standard deviation of residuals from full model fit. 23 Added global variable calc_tstats. 24 Also, added screen display of p-values. 25 Date: 10 May 2000 26 27 */ 28 29 /*---------------------------------------------------------------------------*/ 30 31 #ifndef _NLFIT_HEADER_ 32 #define _NLFIT_HEADER_ 33 34 #include "NLfit_model.h" 35 36 37 /*---------------- global data -------------------*/ 38 static char lbuf[4096] ; 39 static char sbuf[256] ; 40 41 static int calc_tstats = 0; /* set to 1 for calculating t-stats */ 42 43 44 /*---------------------------------------------------------------------------*/ 45 /* 46 Routine to print error message and stop. 47 */ 48 49 void NLfit_error 50 ( 51 char * message /* message to be displayed */ 52 ); 53 54 55 /*---------------------------------------------------------------------------*/ 56 /* 57 Routine to initialize the signal model by defining the number of parameters 58 in the signal model, the default values for the minimum and maximum 59 parameter constraints, and setting the pointer to the function which 60 implements the signal model. 61 */ 62 63 void initialize_signal_model 64 ( 65 NLFIT_MODEL_array * model_array, /* the array of SO models */ 66 char * sname, /* name of signal model (user input) */ 67 vfp * smodel, /* pointer to signal model */ 68 int * p, /* number of parameters in the signal model */ 69 char ** spname, /* signal parameter names */ 70 float * min_sconstr, /* minimum parameter constraints for signal model */ 71 float * max_sconstr /* maximum parameter constraints for signal model */ 72 ); 73 74 75 /*---------------------------------------------------------------------------*/ 76 /* 77 Routine to initialize the noise0 model by defining the number of parameters 78 in the noise model, the default values for the minimum and maximum 79 parameter constraints, and setting the pointer to the function which 80 implements the noise model. 81 */ 82 83 void initialize_noise_model 84 ( 85 NLFIT_MODEL_array * model_array, /* the array of SO models */ 86 char * nname, /* name of noise model (user input) */ 87 vfp * nmodel, /* pointer to noise model */ 88 int * r, /* number of parameters in the noise model */ 89 char ** npname, /* noise parameter names */ 90 float * min_nconstr, /* minimum parameter constraints for noise model */ 91 float * max_nconstr /* maximum parameter constraints for noise model */ 92 ); 93 94 95 /*---------------------------------------------------------------------------*/ 96 /* 97 Routine to calculate the least squares linear regression. 98 */ 99 100 void calc_linear_regression 101 ( 102 matrix x, /* matrix of constants */ 103 vector y, /* observation vector */ 104 vector * b, /* estimated regression coefficients */ 105 float * sse /* error sum of squares */ 106 ); 107 108 109 /*---------------------------------------------------------------------------*/ 110 111 void calc_reduced_model 112 ( 113 int n, /* number of observations */ 114 int r, /* number of parameters in reduced model */ 115 float ** x_array, /* data matrix */ 116 float * y_array, /* observed time series */ 117 float * b_array, /* estimated parameters for reduced model */ 118 float * sse /* error sum of squares */ 119 ); 120 121 122 /*---------------------------------------------------------------------------*/ 123 /* 124 Routine to allocate memory required for calculation of the full model. 125 */ 126 127 void initialize_full_model 128 ( 129 int dimension, /* number of parameters in full model */ 130 int nbest, /* keep nbest vectors from random search */ 131 float *** parameters, /* array of parameter vectors */ 132 float ** sse /* array of error sum of squares */ 133 ); 134 135 136 /*---------------------------------------------------------------------------*/ 137 /* 138 Routine to determine if the parameter vector violates the constraints. 139 */ 140 141 int calc_constraints 142 ( 143 int r, /* number of parameters in the noise model */ 144 int p, /* number of parameters in the signal model */ 145 int nabs, /* use absolute constraints for noise parameters */ 146 float * b_array, /* estimated parameters for the reduced model */ 147 float * min_nconstr, /* minimum parameter constraints for noise model */ 148 float * max_nconstr, /* maximum parameter constraints for noise model */ 149 float * min_sconstr, /* minimum parameter constraints for signal model */ 150 float * max_sconstr, /* maximum parameter constraints for signal model */ 151 float * vertex /* test parameter vector */ 152 ); 153 154 155 /*---------------------------------------------------------------------------*/ 156 /* 157 Calculate the estimated time series using the full model. 158 */ 159 160 void full_model 161 ( 162 vfp nmodel, /* pointer to noise model */ 163 vfp smodel, /* pointer to signal model */ 164 float * gn, /* parameters for noise model */ 165 float * gs, /* parameters for signal model */ 166 int ts_length, /* length of time series data */ 167 float ** x_array, /* independent variable matrix */ 168 float * yhat_array /* output estimated time series */ 169 ); 170 171 172 /*---------------------------------------------------------------------------*/ 173 /* 174 This routine calculates the error sum of squares corresponding to 175 the given vector of parameters for the full model. 176 */ 177 178 float calc_sse 179 ( 180 vfp nmodel, /* pointer to noise model */ 181 vfp smodel, /* pointer to signal model */ 182 int r, /* number of parameters in the noise model */ 183 int p, /* number of parameters in the signal model */ 184 int nabs, /* use absolute constraints for noise parameters */ 185 float * min_nconstr, /* minimum parameter constraints for noise model */ 186 float * max_nconstr, /* maximum parameter constraints for noise model */ 187 float * min_sconstr, /* minimum parameter constraints for signal model */ 188 float * max_sconstr, /* maximum parameter constraints for signal model */ 189 float * b_array, /* estimated parameters for the reduced model */ 190 float * vertex, /* test parameter vector */ 191 int ts_length, /* length of time series array */ 192 float ** x_array, /* independent variable matrix */ 193 float * ts_array /* observed time series */ 194 ); 195 196 197 /*---------------------------------------------------------------------------*/ 198 /* 199 Select and evaluate randomly chosen vectors in the parameter space. 200 */ 201 202 void random_search 203 ( 204 vfp nmodel, /* pointer to noise model */ 205 vfp smodel, /* pointer to signal model */ 206 int r, /* number of parameters in the noise model */ 207 int p, /* number of parameters in the signal model */ 208 int nabs, /* use absolute constraints for noise parameters */ 209 float * min_nconstr, /* minimum parameter constraints for noise model */ 210 float * max_nconstr, /* maximum parameter constraints for noise model */ 211 float * min_sconstr, /* minimum parameter constraints for signal model */ 212 float * max_sconstr, /* maximum parameter constraints for signal model */ 213 int ts_length, /* length of time series array */ 214 float ** x_array, /* independent variable matrix */ 215 float * ts_array, /* observed time series */ 216 float * par_rdcd, /* estimated parameters for the reduced model */ 217 int nrand, /* number of random vectors to generate */ 218 int nbest, /* number of random vectors to keep */ 219 float ** parameters, /* array of best random vectors */ 220 float * response /* array of best sse's */ 221 ); 222 223 224 /*---------------------------------------------------------------------------*/ 225 /* 226 Estimate the parameters for the full model. 227 */ 228 229 void calc_full_model 230 ( 231 vfp nmodel, /* pointer to noise model */ 232 vfp smodel, /* pointer to signal model */ 233 int r, /* number of parameters in the noise model */ 234 int p, /* number of parameters in the signal model */ 235 float * min_nconstr, /* minimum parameter constraints for noise model */ 236 float * max_nconstr, /* maximum parameter constraints for noise model */ 237 float * min_sconstr, /* minimum parameter constraints for signal model */ 238 float * max_sconstr, /* maximum parameter constraints for signal model */ 239 int ts_length, /* length of time series array */ 240 float ** x_array, /* independent variable matrix */ 241 float * ts_array, /* observed time series */ 242 float * par_rdcd, /* estimated parameters for the reduced model */ 243 float sse_rdcd, /* error sum of squares for the reduced model */ 244 int nabs, /* use absolute constraints for noise parameters */ 245 int nrand, /* number of random vectors to generate */ 246 int nbest, /* number of random vectors to keep */ 247 float rms_min, /* minimum rms error required to reject rdcd model */ 248 float * par_full, /* estimated parameters for the full model */ 249 float * sse_full, /* error sum of squares for the full model */ 250 int * novar /* flag for insufficient variation in data */ 251 ); 252 253 254 /*---------------------------------------------------------------------------*/ 255 /* 256 Routine to calculate the partial derivative matrix, evaluated at the 257 estimated parameter location. 258 */ 259 260 void calc_partial_derivatives 261 ( 262 vfp nmodel, /* pointer to noise model */ 263 vfp smodel, /* pointer to signal model */ 264 int r, /* number of parameters in the noise model */ 265 int p, /* number of parameters in the signal model */ 266 float * min_nconstr, /* minimum parameter constraints for noise model */ 267 float * max_nconstr, /* maximum parameter constraints for noise model */ 268 float * min_sconstr, /* minimum parameter constraints for signal model */ 269 float * max_sconstr, /* maximum parameter constraints for signal model */ 270 int ts_length, /* length of time series data */ 271 float ** x_array, /* independent variable matrix */ 272 float * par_full, /* estimated parameters for the full model */ 273 matrix d /* matrix of estimated partial derivatives */ 274 ); 275 276 277 /*---------------------------------------------------------------------------*/ 278 279 void analyze_results 280 ( 281 vfp nmodel, /* pointer to noise model */ 282 vfp smodel, /* pointer to signal model */ 283 int r, /* number of parameters in the noise model */ 284 int p, /* number of parameters in the signal model */ 285 int novar, /* flag for insufficient variation in the data */ 286 float * min_nconstr, /* minimum parameter constraints for noise model */ 287 float * max_nconstr, /* maximum parameter constraints for noise model */ 288 float * min_sconstr, /* minimum parameter constraints for signal model */ 289 float * max_sconstr, /* maximum parameter constraints for signal model */ 290 int ts_length, /* length of time series data */ 291 float ** x_array, /* independent variable matrix */ 292 float * par_rdcd, /* estimated parameters for the reduced model */ 293 float sse_rdcd, /* error sum of squares for the reduced model */ 294 float * par_full, /* estimated parameters for the full model */ 295 float sse_full, /* error sum of squares for the full model */ 296 float * rmsreg, /* root-mean-square for the full regression model */ 297 float * freg, /* f-statistic for the full regression model */ 298 float * rsqr, /* R^2 (coef. of multiple determination) */ 299 float * smax, /* signed maximum of signal */ 300 float * tmax, /* epoch of signed maximum of signal */ 301 float * pmax, /* percentage change due to signal */ 302 float * area, /* area between signal and baseline */ 303 float * parea, /* percentage area between signal and baseline */ 304 float * tpar_full /* t-statistic of the parameters in the full model */ 305 ); 306 307 308 /*---------------------------------------------------------------------------*/ 309 /* 310 Report results for this voxel. 311 */ 312 313 void report_results 314 ( 315 char * nname, /* name of noise model */ 316 char * sname, /* name of signal model */ 317 int r, /* number of parameters in the noise model */ 318 int p, /* number of parameters in the signal model */ 319 char ** npname, /* noise parameter names */ 320 char ** spname, /* signal parameter names */ 321 int ts_length, /* length of time series array */ 322 float * par_rdcd, /* estimated parameters for the reduced model */ 323 float sse_rdcd, /* error sum of squares for the reduced model */ 324 float * par_full, /* estimated parameters for the full model */ 325 float sse_full, /* error sum of squares for the full model */ 326 float * tpar_full, /* t-statistic of the parameters in full model */ 327 float rmsreg, /* root-mean-square for the full regression model */ 328 float freg, /* f-statistic for the full regression model */ 329 float rsqr, /* R^2 (coef. of multiple determination) */ 330 float smax, /* signed maximum of signal */ 331 float tmax, /* epoch of signed maximum of signal */ 332 float pmax, /* percentage change due to signal */ 333 float area, /* area between signal and baseline */ 334 float parea, /* percentage area between signal and baseline */ 335 char ** label /* label output for this voxel */ 336 ); 337 338 339 /*---------------------------------------------------------------------------*/ 340 341 342 #endif /* _NLFIT__HEADER_ */ 343 344