TPCCLIB
|
Header file for libtpcmodel. More...
#include "tpcclibConfig.h"
#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <string.h>
#include <strings.h>
#include <ctype.h>
#include <math.h>
#include <float.h>
#include <sys/time.h>
#include <time.h>
#include "libtpcmisc.h"
Go to the source code of this file.
Data Structures | |
struct | MERTWI |
struct | bobyqa_data |
struct | TGO_POINT |
Macros | |
#define | MAX_PARAMETERS 50 |
#define | MAX_PARAMS MAX_PARAMETERS |
#define | TPCCLIB_MERTWI_NN 312 |
#define | TPCCLIB_MERTWI_A UINT64_C(0xB5026F5AA96619E9) |
#define | MTGA_BEST_MIN_NR 5 |
#define | DEFAULT_SAO2 0.97 |
#define | DEFAULT_P50HB 3.6 |
#define | DEFAULT_P50MB 0.319 |
#define | DEFAULT_NHB 2.7 |
#define | DEFAULT_CHB 150.0 |
#define | DEFAULT_CMB 4.7 |
Typedefs | |
typedef double(* | bobyqa_func) (int n, const double *x, void *func_data) |
Enumerations | |
enum | bobyqa_result { BOBYQA_INVALID_ARGS = -1 , BOBYQA_OUT_OF_MEMORY = -2 , BOBYQA_ROUNDOFF_LIMITED = -3 , BOBYQA_FAIL = -4 , BOBYQA_SUCCESS = 0 , BOBYQA_MINF_MAX_REACHED = 1 , BOBYQA_FTOL_REACHED = 2 , BOBYQA_XTOL_REACHED = 3 , BOBYQA_MAXEVAL_REACHED = 4 , BOBYQA_RELFTOL_REACHED = 5 , BOBYQA_ABSFTOL_REACHED = 6 } |
enum | linefit_method { PEARSON , PERP , LLSQWT , MEDIAN } |
enum | linefit_range { PRESET , EXCLUDE_BEGIN , EXCLUDE_END } |
Functions | |
uint32_t | mertwiSeed32 (void) |
Make uint32_t seed for pseudorandom number generators. | |
uint64_t | mertwiSeed64 (void) |
Make uint64_t seed for pseudorandom number generators. | |
void | mertwiInit (MERTWI *mt) |
void | mertwiInitWithSeed64 (MERTWI *mt, uint64_t seed) |
Initialize the state vector mt[] inside data struct for Mersenne Twister MT19937 pseudorandom number generator using given seed. | |
void | mertwiInitByArray64 (MERTWI *mt, uint64_t init_key[], uint64_t key_length) |
Initialize the state vector mt[] inside data struct for Mersenne Twister MT19937 pseudorandom number generator using given array. | |
uint64_t | mertwiRandomInt64 (MERTWI *mt) |
Generate a random number on [0, 2^64-1]-interval using Mersenne Twister MT19937. | |
int64_t | mertwiRandomInt63 (MERTWI *mt) |
Generate a random number on [0, 2^63-1]-interval using Mersenne Twister MT19937. | |
double | mertwiRandomDouble1 (MERTWI *mt) |
Generate a 64-bit double precision floating point pseudorandom number in the range of [0,1] with uniform distribution using Mersenne Twister MT19937. | |
double | mertwiRandomDouble2 (MERTWI *mt) |
Generate a 64-bit double precision floating point pseudorandom number in the range of [0,1) with uniform distribution using Mersenne Twister MT19937. | |
double | mertwiRandomDouble3 (MERTWI *mt) |
Generate a 64-bit double precision floating point pseudorandom number in the range of (0,1) with uniform distribution using Mersenne Twister MT19937. | |
double | aicSS (double ss, const int n, const int k) |
int | parFreeNr (const int n, double *pLower, double *pUpper) |
Calculate the number of free parameters. | |
int | aicWeights (double *aic, double *w, int n) |
double | aicWeightedAvg (double *w, double *p, int n) |
double | aicModel (double *w, int n) |
bobyqa_result | bobyqb (bobyqa_data *bdata) |
bobyqa_result | bobyqa (int n, int npt, double *x, const double *xl, const double *xu, const double *dx, const double rhoend, double xtol_rel, double minf_max, double ftol_rel, double ftol_abs, int maxeval, int *nevals, double *minf, double(*f)(int n, double *x, void *objf_data), void *objf_data, double *working_space, int verbose) |
int | bobyqa_minimize_single_parameter (bobyqa_data *bdata) |
char * | bobyqa_rc (bobyqa_result rc) |
int | fixed_params (int n, const double *lower, const double *upper, const double *delta) |
int | bobyqa_working_memory_size (int n, int fitted_n, int npt, bobyqa_data *bdata) |
bobyqa_result | bobyqa_set_memory (int n, int fitted_n, int npt, bobyqa_data *bdata, double *wm) |
bobyqa_result | bobyqa_free_memory (bobyqa_data *bdata) |
bobyqa_result | bobyqa_reset_memory (bobyqa_data *bdata) |
void | bobyqa_print (bobyqa_data *bdata, int sw, FILE *fp) |
double | bobyqa_x_funcval (bobyqa_data *bdata, double *x) |
void | bobyqa_xfull (bobyqa_data *bdata) |
bobyqa_result | bobyqa_set_optimization (int full_n, double *x, const double *dx, const double *xl, const double *xu, const double rhoend, double xtol_rel, double minf_max, double ftol_rel, double ftol_abs, int maxeval, double(*f)(int n, double *x, void *objf_data), void *objf_data, int verbose, bobyqa_data *bdata) |
int | bootstrap (int iterNr, double *cLim1, double *cLim2, double *SD, double *parameter, double *lowlim, double *uplim, int frameNr, double *origTac, double *fitTac, double *bsTAC, int parNr, double *weight, double(*objf)(int, double *, void *), char *status, int verbose) |
int | bvls (int key, const int m, const int n, double *a, double *b, double *bl, double *bu, double *x, double *w, double *act, double *zz, int *istate, int *iter, int verbose) |
Bounded-value least-squares method to solve the linear problem A x ~ b , subject to limit1 <= x <= limit2. | |
int | llsqWght (int N, int M, double **A, double *a, double *b, double *weight) |
int | llsqWghtSquared (int N, int M, double **A, double *a, double *b, double *weight) |
int | modelCheckParameters (int par_nr, double *lower_p, double *upper_p, double *test_p, double *accept_p, double *penalty) |
int | modelCheckLimits (int par_nr, double *lower_p, double *upper_p, double *test_p) |
int | fitExpDecayNNLS (double *x, double *y, int n, double fittime, double kmin, double kmax, int pnr, double *a, double *k, int *fnr, int verbose) |
Estimate initial values for sum of exponentials to be fitted on decaying x,y-data. | |
unsigned int | drandSeed (short int seed) |
Make and optionally set the seed for rand(), drand, drandRange, and drandGaussian(). | |
double | gaussdev () |
double | gaussdev2 () |
void | init_gaussdev () |
double | drand () |
int | rand_range (int nr, double *d, double low, double up, int type) |
double | householder_transform (double *v, int N) |
int | householder_hm (double tau, double *vector, double **matrix, int rowNr, int columnNr) |
int | householder_hv (double tau, int size, double *v, double *w) |
double | householder_norm (double *v, int size) |
int | interpolate (double *x, double *y, int nr, double *newx, double *newy, double *newyi, double *newyii, int newnr) |
Linear interpolation and integration. | |
int | finterpolate (float *x, float *y, int nr, float *newx, float *newy, float *newyi, float *newyii, int newnr) |
float version of interpolate(). | |
int | integrate (double *x, double *y, int nr, double *yi) |
int | fintegrate (float *x, float *y, int nr, float *yi) |
float version of integrate(). | |
int | petintegrate (double *x1, double *x2, double *y, int nr, double *newyi, double *newyii) |
int | fpetintegrate (float *x1, float *x2, float *y, int nr, float *newyi, float *newyii) |
Calculates integrals of PET data at frame end times. Float version of petintegrate(). | |
int | interpolate4pet (double *x, double *y, int nr, double *newx1, double *newx2, double *newy, double *newyi, double *newyii, int newnr) |
Interpolate and integrate TAC to PET frames. | |
int | finterpolate4pet (float *x, float *y, int nr, float *newx1, float *newx2, float *newy, float *newyi, float *newyii, int newnr) |
Interpolate and integrate TAC to PET frames. Float version of interpolate4pet(). | |
int | petintegral (double *x1, double *x2, double *y, int nr, double *ie, double *iie) |
Integrate PET TAC data to frame mid times. | |
int | fpetintegral (float *x1, float *x2, float *y, int nr, float *ie, float *iie) |
Integrate PET TAC data to frame mid times. Float version of petintegral(). | |
int | petintegrate2fe (double *x1, double *x2, double *y, int nr, double *e, double *ie, double *iie) |
Integrate PET TAC data to frame end times. | |
int | fpetintegrate2fe (float *x1, float *x2, float *y, int nr, float *e, float *ie, float *iie) |
Integrate PET TAC data to frame end times. Float version of petintegrate2fe(). | |
int | llsqwt (double *x, double *y, int n, double *wx, double *wy, double tol, double *w, double *ic, double *slope, double *nwss, double *sic, double *sslope, double *cx, double *cy) |
int | best_llsqwt (double *x, double *y, double *wx, double *wy, int nr, int min_nr, int mode, double *slope, double *ic, double *nwss, double *sslope, double *sic, double *cx, double *cy, int *bnr) |
int | llsqperp (double *x, double *y, int nr, double *slope, double *ic, double *ssd) |
int | llsqperp3 (double *x, double *y, int nr, double *slope, double *ic, double *ssd) |
int | quadratic (double a, double b, double c, double *m1, double *m2) |
int | medianline (double *x, double *y, int nr, double *slope, double *ic) |
double | least_median_of_squares (double *data, int n) |
Fit a constant (horisontal straight line) to the data by minimising the median of squared residuals. | |
int | least_trimmed_square (double data[], long int n, double *mean, double *variance) |
double | d_kth_smallest (double *data, int n, int k) |
double | dmedian (double *data, int n) |
double | dmean (double *data, int n, double *sd) |
double | dmean_nan (double *data, int n, double *sd, int *vn) |
double | mEstim (double *data, int nr, int iterNr, double cutoff) |
double | huber (double x, double b) |
int | patlak_data (int data_nr, double *i, double *ii, double *c, double *x, double *y) |
int | logan_data (int data_nr, double *i, double *ii, double *c, double *ci, double k2, double *x, double *y) |
int | mtga_best_perp (double *x, double *y, int nr, double *slope, double *ic, double *ssd, int *fnr) |
int | nnls (double **a, int m, int n, double *b, double *x, double *rnorm, double *w, double *zz, int *index) |
int | nnlsWght (int N, int M, double **A, double *b, double *weight) |
int | nnlsWghtSquared (int N, int M, double **A, double *b, double *sweight) |
double | ndtr (double a) |
double | normal_pvalue_2 (double x) |
double | normal_pvalue_1 (double x) |
int | pearson (double *x, double *y, int nr, double *k, double *kSD, double *b, double *bSD, double *r, double *ySD) |
int | pearson2 (double *x, double *y, char *is, int nr, double *k, double *kSD, double *b, double *bSD, double *r, double *ySD) |
int | pearson3 (double *x, double *y, int nr, double *k, double *kSD, double *b, double *bSD, double *r, double *ySD) |
int | pearson4 (double *x, double *y, int nr, double start, double end, double *k, double *kSD, double *b, double *bSD, double *r, double *ySD) |
Calculate slope and intercept of a line and Pearson's correlation coefficient. | |
int | best_pearson (double *x, double *y, int nr, int min_nr, int *first, int *last, double *k, double *kSD, double *b, double *bSD, double *r, double *ySD) |
int | mean (double *x, double *y, int nr, double *xmean, double *xsd, double *ymean, double *ysd) |
int | regr_line (double *x, double *y, int n, double *m, double *c) |
int | highest_slope (double *x, double *y, int n, int slope_n, double *m, double *c, double *xi, double *xh) |
int | highest_slope_after (double *x, double *y, int n, int slope_n, double x_start, double *m, double *c, double *xi, double *xh) |
int | powell (double *p, double *delta, int parNr, double ftol, int *iterNr, double *fret, double(*_fun)(int, double *, void *), void *fundata, int verbose) |
int | qrLSQ (double **mat, double *rhs, double *sol, const unsigned int rows, const unsigned int cols, double *r2) |
QR least-squares solving routine. | |
int | qr (double **A, int m, int n, double *B, double *X, double *rnorm, double *tau, double *res, double **wws, double *ws) |
int | qr_decomp (double **a, int M, int N, double *tau, double **cchain, double *chain) |
int | qr_solve (double **QR, int M, int N, double *tau, double *b, double *x, double *residual, double *resNorm, double **cchain, double *chain) |
int | qr_weight (int N, int M, double **A, double *b, double *weight, double *ws) |
int | qrLH (const unsigned int m, const unsigned int n, double *a, double *b, double *x, double *r2) |
Solve over-determined least-squares problem A x ~ b using successive Householder rotations. | |
int | runs_test (double *data1, double *data2, int N, double alpha, double *p) |
int | residuals (double *d1, double *d2, int Nr, int *Rnr, int *Nminus, int *Nplus) |
int | mrl_between_tacs (double y1[], double y2[], int n) |
void | random_shuffle (int *array, int n) |
void | randperm (int *array, int n, int a) |
double | simplex (double(*_fun)(double *), int parNr, double *par, double *delta, double maxerr, int maxiter, int verbose) |
int | simC3s (double *t, double *ca, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double *ct, double *cta, double *ctb, double *ctc) |
int | simC3p (double *t, double *ca, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double *ct, double *cta, double *ctb, double *ctc) |
int | simC3vs (double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb) |
int | simC3vp (double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb) |
int | simC2l (double *t, double *ca, int nr, double k1, double k2, double k3, double kLoss, double *ct, double *cta, double *ctb) |
int | simC2vl (double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double kL, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctab, double *ctvb) |
int | simC3vpKLoss (double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double kLoss, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb) |
int | simRTCM (double *t, double *cr, int nr, double R1, double k2, double k3, double k4, double *ct, double *cta, double *ctb) |
int | simSRTM (double *t, double *cr, int nr, double R1, double k2, double BP, double *ct) |
int | simTRTM (double *t, double *cr, int nr, double R1, double k2, double k3, double *ct) |
int | simHuangmet (double *t, double *ctot, int nr, double k01, double k12, double k21, double k03, double k34, double k43, double *c0, double *c1, double *c3) |
int | simTPCMOD0009c (double *t, double *ctot, int nr, double km, double k1m, double k2m, double k3m, double k4m, double *ca, double *cm) |
int | simMBF (double *t, double *ci, int nr, double k1, double k2, double Vfit, double *ct) |
int | simC1 (double *t, double *ca, int nr, double k1, double k2, double *ct) |
int | simC3DIvs (double *t, double *ca1, double *ca2, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double k1b, double k2b, double f, double vb, double fa, double *scpet, double *sct1, double *sct2, double *sct3, double *sct1b, double *sctab, double *sctvb) |
int | simC4DIvp (double *t, double *ca1, double *ca2, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double k7, double km, double k1b, double k2b, double f, double vb, double fa, double *scpet, double *sct1, double *sct2, double *sct3, double *sct1b, double *sctab, double *sctvb, int verbose) |
int | simC4DIvs (double *t, double *ca1, double *ca2, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double k7, double km, double k1b, double k2b, double f, double vb, double fa, double *scpet, double *sct1, double *sct2, double *sct3, double *sct1b, double *sctab, double *sctvb, int verbose) |
int | simDispersion (double *x, double *y, int n, double tau1, double tau2, double *tmp) |
int | simOxygen (double *t, double *ca1, double *ca2, double *ca1i, double *ca2i, const int n, const double k1a, const double k2a, const double km, const double k1b, const double k2b, const double vb, const double fa, double *scpet, double *sct1, double *sct2, double *sctab, double *sctvb1, double *sctvb2, double *scvb1, double *scvb2, const int verbose) |
double | mo2k1k2 (const double OER, const double SaO2, const double p50Hb, const double p50Mb, const double nHb, const double cHb, const double cMb, const int verbose) |
Calculates K1/k2 ratio for [O-15]O2 in muscle, based on OER. | |
double | mo2pO2 (const double OER, const double K1k2, const double SaO2, const double p50Mb, const double cHb, const double cMb, const int verbose) |
int | tgo (double *lowlim, double *uplim, double(*objf)(int, double *, void *), void *objfData, int dim, int neighNr, double *fmin, double *gmin, int samNr, int tgoNr, int verbose) |
void | tgoRandomParameters (TGO_POINT *p, int parNr, int sNr, double *low, double *up) |
void | tgoRandomParametersST (TGO_POINT *p, int parNr, int sNr, double *low, double *up) |
int | nlopt1D (double(*_fun)(double, void *), void *_fundata, double x, double xl, double xu, double delta, double tol, const int maxeval, double *nx, double *nf, int verbose) |
Variables | |
long int | GAUSSDEV_SEED |
int | TGO_SQUARED_TRANSF |
int | TGO_LOCAL_INSIDE |
int | TGO_LOCAL_OPT |
Header file for libtpcmodel.
Definition in file libtpcmodel.h.
#define DEFAULT_CHB 150.0 |
Hemoglobin concentration in blood (mg/g)
Definition at line 962 of file libtpcmodel.h.
#define DEFAULT_CMB 4.7 |
Myoglobin concentration in muscle (mg/g)
Definition at line 964 of file libtpcmodel.h.
#define DEFAULT_NHB 2.7 |
Hill coefficient n for hemoglobin
Definition at line 960 of file libtpcmodel.h.
#define DEFAULT_P50HB 3.6 |
Half-saturation pressure p50 (kPa) for hemoglobin
Definition at line 956 of file libtpcmodel.h.
#define DEFAULT_P50MB 0.319 |
Half-saturation pressure p50 (kPa) p50 for myoglobin
Definition at line 958 of file libtpcmodel.h.
#define DEFAULT_SAO2 0.97 |
Arterial oxygen saturation fraction
Definition at line 954 of file libtpcmodel.h.
#define MAX_PARAMETERS 50 |
#define MAX_PARAMS MAX_PARAMETERS |
Max nr of parameters
Definition at line 35 of file libtpcmodel.h.
#define MTGA_BEST_MIN_NR 5 |
Min nr of points in MTGA line fit
Definition at line 664 of file libtpcmodel.h.
Referenced by img_logan(), img_patlak(), and mtga_best_perp().
#define TPCCLIB_MERTWI_A UINT64_C(0xB5026F5AA96619E9) |
Mersenne Twister required constant
Definition at line 43 of file libtpcmodel.h.
Referenced by mertwiInit(), and mertwiRandomInt64().
#define TPCCLIB_MERTWI_NN 312 |
Mersenne Twister state vector length
Definition at line 41 of file libtpcmodel.h.
Referenced by mertwiInit().
typedef double(* bobyqa_func) (int n, const double *x, void *func_data) |
BOBYQA objective function
Definition at line 92 of file libtpcmodel.h.
enum bobyqa_result |
BOBYQA return codes
Definition at line 77 of file libtpcmodel.h.
enum linefit_method |
MTGA line fit method
Definition at line 667 of file libtpcmodel.h.
enum linefit_range |
MTGA line fit range
Definition at line 681 of file libtpcmodel.h.
|
extern |
Calculates a value describing the relative goodness of models, based on an array of model weights.
w | Array of weights |
n | Length of array |
Definition at line 132 of file aic.c.
|
extern |
Computation of AICc in the special case of sum-of-squares optimization from the SS, nr of fitted points and nr of fitted parameters.
If variance is different between the data points, weighted SS must be given.
ss | Sum-of-Squares of the fit |
n | Sample size, i.e. nr of fitted data points |
k | Number of fitted model parameters; AICc calculation is valid only when (n-k)>1. |
Definition at line 20 of file aic.c.
Referenced by imgNoiseTemplate().
|
extern |
Computation of the Akaike weighted model parameter average. Requires arrays of AIC weight values, and corresponding parameter values.
w | Array of weights |
p | Array of parameters |
n | Lengths of arrays |
Definition at line 108 of file aic.c.
|
extern |
Computation of the Akaike weights for model averaging. Requires an array of AIC values, and an output array for weights.
aic | Array of AICs |
w | Array of weights (output) |
n | Lengths of arrays |
Definition at line 74 of file aic.c.
|
extern |
Finds the best least-squares line to (x,y)-data, leaving points out either from the beginning (mode=0) or from the end (mode=1).
x | Plot x axis values. |
y | Plot y axis values. |
wx | Weighting factors for x. |
wy | Weighting factors for y. |
nr | Nr of plot data points. |
min_nr | Min nr of points to use in the fit; must be >=4. |
mode | Leave out points from beginning (0) or from end (1). |
slope | Slope is returned in here. |
ic | Y axis intercept is returned in here. |
nwss | sqrt(WSS)/wsum is returned here. |
sslope | Expected sd of slope at calculated points. |
sic | Expected sd of intercept at calculated points. |
cx | Calculated x data points. |
cy | Calculated y data points. |
bnr | Number of points in the best fit (incl data with w=0). |
Definition at line 294 of file llsqwt.c.
|
extern |
Find the best linear fit to double data (x[], y[]) with nr points.
Data may contain NaN's, which are not used.
x | Data x values |
y | Data y values |
nr | Number of data sample values |
min_nr | Minimum nr of data points to use |
first | Index [0..last-2] of the first point to use initially, and after fitting |
last | Index [first+1..nr-1] of the last point to use initially, and after fitting |
k | slope |
kSD | S.D. of slope |
b | y axis intercept |
bSD | S.D. of y axis intercept |
r | Pearson's correlation coefficient, r |
ySD | Residual variance of y values |
Definition at line 252 of file pearson.c.
|
extern |
BOBYQA is a derivative-free optimization code with constraints. This function interface is close to the one published by Powell.
n | N must be set to the number of variables (and length of arrays below) |
npt | NPT is the number of interpolation conditions. Its value must be in the interval [N+2,(N+1)(N+2)/2]. Choices that exceed 2*N+1 are not recommended. If set to zero, then bobyqa uses 2*N+1 as default. However, if very high precision is required, relatively high npt seems to be needed. |
x | Initial values of the variables must be set in X(1),X(2),...,X(N). They will be changed to the values that give the least calculated F. |
xl | Lower bounds for X[]: For I=1,2,...,N, XL(I) and XU(I) must provide the lower and upper bounds, respectively, on X(I). The construction of quadratic models requires XL(I) to be strictly less than XU(I) for each I. Further, the contribution to a model from changes to the I-th variable is damaged severely by rounding errors if XU(I)-XL(I) is too small. If xl[i]==xu[i]==x[i], then parameter x[i] is fixed to x[i] and left out of actual bobyqa procedure; parameter(s) can also be fixed by setting dx[i]=0. |
xu | Upper bounds for X[]: For I=1,2,...,N, XL(I) and XU(I) must provide the lower and upper bounds, respectively, on X(I). The construction of quadratic models requires XL(I) to be strictly less than XU(I) for each I. Further, the contribution to a model from changes to the I-th variable is damaged severely by rounding errors if XU(I)-XL(I) is too small. If xl[i]==xu[i]==x[i], then parameter x[i] is fixed to x[i] and left out of actual bobyqa procedure; parameter(s) can also be fixed by setting dx[i]=0. |
dx | Initial step sizes. dx should be set to 1/10 of probable distance from initial guess where the final solution could be; set dx to zero for those parameters that are to be fixed to their initial value. Too high values may lead to error. |
rhoend | Rhoend must be set to the final value of a trust region radius, as a (positive) fraction of dx values given above, for example 1.0E-06. Smaller value adds precision and increases calculation time. Enter <=0 to use default value, calculated using xtol_rel. |
xtol_rel | Relative X tolerance |
minf_max | Stopping rule: maximal allowed function value; when reached, stops. |
ftol_rel | Stopping rule: relative tolerance to function value |
ftol_abs | Stopping rule: absolute tolerance to function value |
maxeval | Stopping rule: max nr of function evaluations |
nevals | Number of function evaluations; this function will initialize this to zero, and in the end returns the number. Enter NULL if not needed. |
minf | Pointer where minimized function value will be written |
f | Objective function has to be provided by the user. It must set F to the value of the objective function for the current values of the variables X(1),X(2),...,X(N), which are generated automatically in a way that satisfies the bounds given in XL and XU. |
objf_data | Pointer to data that is carried over to the function. It is not needed or used by bobyqa. Enter NULL if not needed. |
working_space | The array W will be used for working space. Its length can be calculated using function bobyqa_working_memory_size(). Enter NULL to leave allocating to this function. |
verbose | Verbose level; if zero, then nothing is printed to stderr or stdout |
Definition at line 40 of file bobyqa.c.
Referenced by tgo().
|
extern |
Free the working memory and memory for data arrays and matrices in BOBYQA structure that was allocated by bobyqa_allocate_memory().
bdata | Pointer to bobyqa struct |
Definition at line 592 of file bobyqa.c.
Referenced by bobyqa(), and bobyqa_set_memory().
|
extern |
Local one-dimensional minimization, meant to be called from bobyqa() which does not otherwise work with 1-D problems.
Definition at line 210 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
For testing: print BOBYQA data struct variables to specified file pointer.
bdata | Pointer to bobyqa structure. |
sw | Switch that tells what to print: 0=All; 1=parameters; 2=array sizes; 3=full parameter list; 4=fitted parameters, limits and deltas. |
fp | File pointer where contents are printed; usually stdout. |
Definition at line 640 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
Description of bobyqa return code as text string.
Definition at line 381 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
Initiate BOBYQA data struct variables before calling bobyqb(). Can be used to reset data variables before calling bobyqb() again with the same problem size, without need to call bobyqa_free_memory() and bobyqa_set_memory() first.
bdata | Pointer to bobyqa struct |
Definition at line 614 of file bobyqa.c.
Referenced by bobyqa_set_memory().
|
extern |
Setup working memory and memory for data arrays and matrices in BOBYQA structure. Double memory is either allocated inside this function or preallocated memory can be given as argument.
n | Number of all model parameters. |
fitted_n | Number of fitted (not fixed) model parameters. |
npt | NPT. |
bdata | Pointer to bobyqa structure. |
wm | Enter pointer to preallocated working memory; NULL to allocate it here. |
Definition at line 501 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
Set optimization parameters and limits in BOBYQA data struct.
full_n | Length of full parameter list. |
x | Pointer to array of full parameters. |
dx | Pointer to array of initial step sizes (deltas). |
xl | Pointer to array of lower limits. |
xu | Pointer to array of upper limits. |
rhoend | Rhoend; enter <=0 to use default value, calculated using xtol_rel. |
xtol_rel | Relative X tolerance. |
minf_max | Stopping rule: maximal allowed function value; when reached, stops. |
ftol_rel | Stopping rule: relative tolerance to function value. |
ftol_abs | Stopping rule: absolute tolerance to function value. |
maxeval | Stopping rule: max nr of function evaluations. |
f | Pointer to object function |
objf_data | Pointer to data that is carried over to the object function. Enter NULL if not needed. |
verbose | Verbose level for all BOBYQA functions; if zero, then nothing is printed to stderr or stdout |
bdata | Pointer to bobyqa structure. |
Definition at line 781 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
Calculate the required size of memory to allocate for bobyqa().
n | Number of all model parameters |
fitted_n | Number of fitted (not fixed) model parameters |
npt | NPT |
bdata | Enter pointer to bobyqa struct to set array sizes, otherwise enter NULL |
Definition at line 441 of file bobyqa.c.
Referenced by bobyqa(), and bobyqa_set_memory().
|
extern |
Calculate the objective function value with current scaled parameters values stored in specified array. This function back-scales the parameters to the full parameter list and calls the objective function with full set of parameters, that is, including also fixed parameters which are saved in bobyqa data struct. Also nevals is incremented.
Definition at line 751 of file bobyqa.c.
Referenced by bobyqa_minimize_single_parameter(), bobyqa_prelim(), bobyqa_rescue(), and bobyqb_calc_with_xnew().
|
extern |
Add the fitted parameters (x[]) to the full parameter list (xfull[]). xfull[] will contain back-scaled parameter values. x[] is not modified.
Definition at line 769 of file bobyqa.c.
Referenced by bobyqa(), and bobyqb_xupdate().
|
extern |
bobyqb (called from bobyqa).
bdata | Data struct for BOBYQA subroutines |
Definition at line 1622 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
Bootstrap method.
Original weights are assumed to be inversely proportional to variance. Square root is used, because bootstrap assumes them to be proportional to standard deviation. If only standard deviation is needed then cLim1 and cLim2 can be set to be NULL, and if only the confidence limits are wanted then SD can be set to be NULL. This function will not set seed for random number generator, therefore, make sure that it is set in your program, for example with srand(time(NULL));.
iterNr | Bootstrap iteration number (>=100), set to zero to use the default (200). |
cLim1 | Vector to put the lower confidence limits to; NULL if not needed. |
cLim2 | Vector to put the upper confidence limits to; NULL if not needed. |
SD | Vector to put the standard deviations to; NULL if not needed. |
parameter | Best parameter estimates (preserved by this function). |
lowlim | Lower limits for the parameters. |
uplim | Upper limits for the parameters. |
frameNr | Nr of samples in tissue TAC data. |
origTac | Measured tissue TAC values (not modified; local copy is used). |
fitTac | Best fitted tissue TAC values (not modified; local copy is used). |
bsTac | Pointer to (empty) tissue TAC vector, where bootstrapped TACs will be written in each bootstrap iteration, and which will be used by objective function to calculate WSS. |
parNr | Nr of parameters. |
weight | sample weights. |
objf | The object function. |
status | Pointer to a string (allocated for at least 64 chars) where error message or other execution status will be written; enter NULL, if not needed. |
verbose | Verbose level; if zero, then nothing is printed to stderr or stdout. |
Definition at line 55 of file bootstrap.c.
|
extern |
Bounded-value least-squares method to solve the linear problem A x ~ b , subject to limit1 <= x <= limit2.
This routine is based on the text and Fortran code in C.L. Lawson and R.J. Hanson, Solving Least Squares Problems, Prentice-Hall, Englewood Cliffs, New Jersey, 1974, and Fortran codes by R.L. Parker and P.B. Stark, and by J. Burkardt.
key | Enter 0 to solve the problem from scratch, or <>0 to initialize the routine using caller's guess about which parameters are active within their bounds, which are at their lower bounds, and which at their upper bounds, as supplied in the array istate ('warm start'). When key <> 0, the routine initially sets the active components to the averages of their upper and lower bounds. |
m | Number of samples in matrix A and the length of vector b. |
n | Number of parameters in matrix A and the length of vector x. The n must not be > m unless at least n-m variables are non-active (set to their bounds).. |
a | Pointer to matrix A; matrix must be given as an n*m array, containing n consecutive m-length vectors. Contents of A are modified in this routine. |
b | Pointer to vector b of length m. Contents of b are modified in this routine. |
bl | Array BL[0..n-1] of lower bounds for parameters. |
bu | Array BU[0..n-1] of upper bounds for parameters. |
x | Pointer to the result vector x of length n. |
w | Pointer to an array of length n, which will be used as working memory, and the minimum 2-norm || a.x-b ||, (R^2), will be written in w[0]. |
act | Pointer to an array of length m*(n+2), or m*(m+2) if m<n, which will be used as working memory. |
zz | Pointer to an array of length m, to be used as working memory. |
istate | Pointer to an integer array of length n+1. If parameter key <>0, then the last position istate[n] must contain the total number of components at their bounds (nbound, the ‘bound variables’). The absolute values of the first nbound entries of istate[] are the indices of these ‘bound’ components of x[]. The sign of istate[0.. nbound-1] entries indicates whether x[|istate[i]|] is at its upper (positive) or lower (negative) bound. Entries istate[nbound..n-1] contain the indices of the active components. |
iter | Number of performed iterations is returned in this variable. Maximum nr of iterations is set with this variable, or set it to 0 to use the default maximum, 3*n. |
verbose | Verbose level; if zero, then nothing is printed to stderr or stdout |
Definition at line 24 of file bvls.c.
|
extern |
Returns the kth smallest value in data[0..n-1]. Array is partially sorted. Algorithm is based on the book Wirth N. Algorithms + data structures = programs. Englewood Cliffs, Prentice-Hall, 1976.
data | Pointer to data; array is partially sorted. |
n | Length of data array. |
k | kth smallest value will be returned. |
Definition at line 15 of file median.c.
Referenced by dmedian().
|
extern |
Returns the mean in array data[0..n-1], and optionally calculates also the (sample) standard deviation of the mean.
data | Pointer to data; data is not changed in any way. |
n | Length of data array. |
sd | Pointer to variable where SD will be written; enter NULL if not needed. |
Definition at line 73 of file median.c.
|
extern |
Returns the mean in array data[0..n-1], and optionally calculates also the standard deviation of the mean. Data may contain missing samples marked as NaNs.
data | Pointer to data; data is not changed in any way. |
n | Length of data array. |
sd | Pointer to variable where SD will be written; enter NULL if not needed. |
vn | Pointer to variable where number of valid (not NaN) samples will be written; enter NULL if not needed. |
Definition at line 104 of file median.c.
|
extern |
Returns the median in array data[0..n-1]. Array is partially sorted. Algorithm is based on the book Wirth N. Algorithms + data structures = programs. Englewood Cliffs, Prentice-Hall, 1976.
data | Pointer to data; array is partially sorted. |
n | Length of data array. |
Definition at line 48 of file median.c.
Referenced by dftRobustMinMaxTAC(), least_trimmed_square(), and mEstim().
|
extern |
Alternative function to rand() which returns a double precision floating point number in the range of [0,1].
Definition at line 142 of file gaussdev.c.
Referenced by gaussdev(), gaussdev2(), rand_range(), tgoRandomParameters(), and tgoRandomParametersST().
|
extern |
Make and optionally set the seed for rand(), drand, drandRange, and drandGaussian().
Uses microseconds from the computer clock and process ID to reduce the chance of getting the same seed for simultaneously executing program threads and instances.
seed | Also sets seed with srand (1) or not (0) |
Definition at line 63 of file gaussdev.c.
Referenced by tgo().
|
extern |
float version of integrate().
Linear integration from time 0 to x[0..nr-1]. If x[0] is >0 and x[0]<=(x[1]-x[0]), then the beginning is interpolated from it to (0,0).
x | Original x values; duplicates are not allowed, all must be >=0. Data must be sorted by ascending x |
y | Original y values |
nr | Nr of values |
yi | Array for integrals |
Definition at line 300 of file integr.c.
|
extern |
float version of interpolate().
It is assumed that both original and new interpolated data represent the actual values at specified time points (not from framed data). The integration is calculated dot-by-dot.
If NULL is specified for newy[], newyi[] and/or newyii[], their values are not calculated.
Original data and new x values must be sorted by ascending x. Subsequent x values can have equal values, enabling the use of step functions. Negative x (time) values can be processed. If necessary, the data is extrapolated assuming that: 1) y[inf]=y[nr-1] and 2) if x[0]>0, y[0]>0 and x[0]<=x[1]-x[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x | Input data x (time) values |
y | Input data y values |
nr | Number of values in input data |
newx | Output data x values |
newy | Interpolated (extrapolated) y values; NULL can be given if not needed |
newyi | Integrals; NULL can be given if not needed |
newyii | 2nd integrals; NULL can be given if not needed |
newnr | Nr of values in output data |
Definition at line 161 of file integr.c.
Referenced by finterpolate4pet(), and imgTimeIntegral().
|
extern |
Interpolate and integrate TAC to PET frames. Float version of interpolate4pet().
It is assumed, that original data is not from framed data, but that the values represent the actual value at specified time point, which allows the integration to be calculated dot-by-dot.
If NULL is specified for *newy, *newyi and/or *newyii, their values are not calculated.
Original data and new x values must be sorted by ascending x. Subsequent x values can have equal values, enabling the use of step functions. PET frames can overlap, but interpolation may then be slower. Negative x (time) values can be processed. If necessary, the data is extrapolated assuming that 1) y[inf]=y[nr-1] and 2) if x[0]>0 and y[0]>0 and x[0]>x[1]-x[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x | Times of original data |
y | Values of original data |
nr | Number of original data values |
newx1 | PET frame start times; frames may overlap |
newx2 | PET frame end times; frames may overlap |
newy | Mean value during PET frame, or NULL if not needed; calculation may be faster if newyi is calculated too |
newyi | Integral at frame mid time, or NULL if not needed |
newyii | 2nd integral at frame mid time, or NULL if not needed |
newnr | Number of PET frames |
Definition at line 646 of file integr.c.
|
extern |
Estimate initial values for sum of exponentials to be fitted on decaying x,y-data.
x | Pointer to array of x values; must not contain NaNs; data is not changed; samples must be sorted by increasing x. |
y | Pointer to array of y values; must not contain NaNs; data is not changed. |
n | Nr of x and y samples (x and y array lengths); data is not changed. |
fittime | Fittime is usually set to <=0 or a very high value to start the search for the best line fit from the last x,y sample towards the first sample; However, to exclude the end phase you may want to set fittime to include only certain time range from the beginning. |
kmin | Minimum eigenvalue, must be >0; for example, 1.0E-06. |
kmax | Maximum eigenvalue; for example, 1.0E+03. |
pnr | Max number of exp functions to save; length of a[] and k[] arrays. |
a | Pointer to an array of length pnr where exp function coefficients will be written; enter NULL if not needed. |
k | Pointer to an array of length pnr where exp function eigenvalues will be written; enter NULL if not needed. |
fnr | The number of fitted exponentials will be written in here; note that this number may be higher than pnr; enter NULL if not needed. |
verbose | Verbose level; if zero, then nothing is printed to stderr or stdout |
Definition at line 100 of file constraints.c.
|
extern |
Determine how many of model parameters are fixed.
n | Total nr of parameters |
lower | List of lower parameter limits |
upper | List of upper parameter limits |
delta | List of parameter deltas; enter NULL if not available |
Definition at line 418 of file bobyqa.c.
Referenced by bobyqa().
|
extern |
Integrate PET TAC data to frame mid times. Float version of petintegral().
Any of output arrays may be set to NULL if that is not needed. Frames must be in ascending time order. Gaps and small overlap are allowed. If x1[0]>0 and x1[0]<=x2[0]-x1[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x1 | frame start times |
x2 | frame end times |
y | avg value during frame |
nr | number of frames |
ie | integrals at frame mid time |
iie | 2nd integrals at frame mid time |
Definition at line 838 of file integr.c.
Referenced by img_logan(), and img_patlak().
|
extern |
Calculates integrals of PET data at frame end times. Float version of petintegrate().
Data does not have to be continuous, but it must be increasing in time. For faster performance in repetitive calls, allocate memory for integral[] even if it is not needed. If x1[0] is >0 AND x1[0] is <=(x2[0]-x1[0]), then the beginning is interpolated from (0, 0) to (x1[0], y1[0]*x1[0]/x) where x is midtime of the first frame.
x1 | Array of frame start times |
x2 | Array of frame end times |
y | Array of y values (avg during each frame) |
nr | Nr of frames |
newyi | Output: integral values at frame end times, or NULL |
newyii | Output: 2nd integral values at frame end times, or NULL |
Definition at line 418 of file integr.c.
|
extern |
Integrate PET TAC data to frame end times. Float version of petintegrate2fe().
Any of output arrays may be set to NULL if that is not needed. Frames must be in ascending time order. Gaps and small overlap are allowed. If x1[0]>0 and x1[0]<=x2[0]-x1[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x1 | frame start times |
x2 | frame end times |
y | avg value during frame |
nr | number of frames |
e | values at frame end time |
ie | integrals at frame end time |
iie | 2nd integrals at frame end time |
Definition at line 974 of file integr.c.
|
extern |
Applies the polar form of Box-Müller transform to produce pseudo-random numbers with Gaussian (normal) distribution which has a zero mean and standard deviation of one. Box GEP, Muller ME. A note on the generation of random normal deviates. Annals of Mathematical Statistics, Volume 29, Issue 2, 1958, 610-611. Available from JSTOR https://www.jstor.org/
Definition at line 30 of file gaussdev.c.
|
extern |
Applies the polar form of Box-Müller transform to produce pseudo-random numbers with Gaussian (normal) distribution which has a zero mean and standard deviation of one. This function does never set seed, like gaussdev() does, therefore set seed for random number generator before first calling this routine, for example with srand(time(NULL));
Box GEP, Muller ME. A note on the generation of random normal deviates. Annals of Mathematical Statistics, Volume 29, Issue 2, 1958, 610-611. Available from JSTOR https://www.jstor.org/
Definition at line 112 of file gaussdev.c.
|
extern |
Finds the regression line with the highest slope for x,y data.
x | An array of x axis values |
y | An array of y axis values |
n | The number of values in x and y arrays |
slope_n | The number of samples used to fit the line |
m | Pointer where calculated slope is written; NULL if not needed |
c | Pointer where calculated y axis intercept is written; NULL if not needed |
xi | Pointer where calculated x axis intercept is written; NULL if not needed |
xh | Pointer where the place (x) of the highest slope is written; NULL if not needed |
Definition at line 424 of file pearson.c.
Referenced by dftFixPeak().
|
extern |
Finds the regression line with the highest slope for x,y data after specified x.
x | An array of x axis values. |
y | An array of y axis values. |
n | The number of values in x and y arrays. |
slope_n | The number of samples used to fit the line. |
x_start | Estimation start x value, samples with smaller x are ignored; can usually be set to zero. |
m | Pointer where calculated slope is written; NULL if not needed. |
c | Pointer where calculated y axis intercept is written; NULL if not needed. |
xi | Pointer where calculated x axis intercept is written; NULL if not needed. |
xh | Pointer where the place (x) of the highest slope is written; NULL if not needed. |
Definition at line 475 of file pearson.c.
|
extern |
Applies a householder transformation defined by vector "vector" and scalar tau to the left-hand side of the matrix. (I - tau vector vector^T)*matrix The result of the transform is stored in matrix.
tau | Coefficient defining householder transform. |
vector | Vector defining householder transform (of size rowNr). |
matrix | the matrix that is to be transformed. |
rowNr | Nr of rows in matrix. |
columnNr | Nr of columns in matrix. |
Definition at line 68 of file hholder.c.
Referenced by qr_decomp().
|
extern |
Applies a householder transformation defined by vector v and coefficient tau to vector w w = (I - tau v v^T) w.
tau | Coefficient defining householder transform. |
size | Size of vectors v and w. |
v | Vector v. |
w | Vector w. |
Definition at line 107 of file hholder.c.
Referenced by qr_solve().
|
extern |
|
extern |
This function prepares a Householder transformation P = I - tau h h^T which can be used to zero all the elements of the input vector except the first one that will get value beta. On output the elements 1 - size-1 of the vector h are stored in locations vector[1] - vector[size-1] of the input vector and value of beta is stored in location vector[0].
v | The N-vector to be transformed. |
N | size of the vector. |
Definition at line 23 of file hholder.c.
Referenced by qr_decomp().
|
extern |
|
extern |
Initiate random number generator for gaussdev()
Definition at line 92 of file gaussdev.c.
Referenced by gaussdev().
|
extern |
Linear integration from time 0 to x[0..nr-1]. If x[0] is >0 and x[0]<=(x[1]-x[0]), then the beginning is interpolated from it to (0,0).
x | Original x values; duplicates are not allowed, all must be >=0. Data must be sorted by ascending x |
y | Original y values |
nr | Nr of values |
yi | Array for integrals |
Definition at line 271 of file integr.c.
Referenced by dftReadReference().
|
extern |
Linear interpolation and integration.
It is assumed that both original and new interpolated data represent the actual values at specified time points (not from framed data). The integration is calculated dot-by-dot.
If NULL is specified for newy[], newyi[] and/or newyii[], their values are not calculated.
Original data and new x values must be sorted by ascending x. Subsequent x values can have equal values, enabling the use of step functions. Negative x (time) values can be processed. If necessary, the data is extrapolated assuming that: 1) y[inf]=y[nr-1] and 2) if x[0]>0, y[0]>0 and x[0]<=x[1]-x[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x | Input data x (time) values |
y | Input data y values |
nr | Number of values in input data |
newx | Output data x values |
newy | Interpolated (extrapolated) y values; NULL can be given if not needed |
newyi | Integrals; NULL can be given if not needed |
newyii | 2nd integrals; NULL can be given if not needed |
newnr | Nr of values in output data |
Definition at line 28 of file integr.c.
Referenced by bfIrr2TCM(), bfRadiowater(), dftDivideFrames(), dftDoubleFrames(), dftInterpolate(), dftInterpolateInto(), dftReadinput(), dftTimeIntegral(), and interpolate4pet().
|
extern |
Interpolate and integrate TAC to PET frames.
It is assumed, that original data is not from framed data, but that the values represent the actual value at specified time point, which allows the integration to be calculated dot-by-dot.
If NULL is specified for *newy, *newyi and/or *newyii, their values are not calculated.
Original data and new x values must be sorted by ascending x. Subsequent x values can have equal values, enabling the use of step functions. PET frames can overlap, but interpolation may then be slower. Negative x (time) values can be processed. If necessary, the data is extrapolated assuming that 1) y[inf]=y[nr-1] and 2) if x[0]>0 and y[0]>0 and x[0]>x[1]-x[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x | Times of original data |
y | Values of original data |
nr | Number of original data values |
newx1 | PET frame start times; frames may overlap |
newx2 | PET frame end times; frames may overlap |
newy | Mean value during PET frame, or NULL if not needed; calculation may be faster if newyi is calculated too |
newyi | Integral at frame mid time, or NULL if not needed |
newyii | 2nd integral at frame mid time, or NULL if not needed |
newnr | Number of PET frames |
Definition at line 510 of file integr.c.
Referenced by bfIrr2TCM(), bfRadiowater(), dftAutointerpolate(), dftInterpolate(), dftInterpolateForIMG(), dftInterpolateInto(), img_k1_using_ki(), img_logan(), and img_patlak().
|
extern |
Fit a constant (horisontal straight line) to the data by minimising the median of squared residuals.
The algorithm is described in P.J. Rousseeuw: Least Median of Squares Regression, Journal of the American Statistical Association, Vol. 79, No. 388 (1984), 871-880.
data | Data array |
n | Number of data values |
Definition at line 21 of file lms.c.
|
extern |
Least trimmed squares estimates for univariate location and variance. Data samples are expected to be truly real valued (i.e too many samples having the same value might lead to problems. The algorithm (exact) is described in P.J. Rousseeuw and A.M. Leroy: Robust Regression and Outlier Detection John Wiley & Sons 1987.
data | Input vector of n sample values; data samples are expected to be truly real valued (i.e too many samples having the same value might lead to problems. |
n | Number of samples |
mean | Output: Mean of sample values |
variance | Output: Variance of sample values |
Definition at line 27 of file lts.c.
|
extern |
Simple non-iterative perpendicular line fitting. This function is fully based on the article [1].
References:
x | Coordinates of data points (dimension nr). |
y | Coordinates of data points (dimension nr). |
nr | Number of data points. |
slope | Estimated slope. |
ic | Estimated intercept. |
ssd | Sum of squared distances / nr. |
Definition at line 382 of file llsqwt.c.
Referenced by img_logan(), img_patlak(), llsqperp3(), and mtga_best_perp().
|
extern |
Simple non-iterative perpendicular line fitting. This version of the function accepts data that contains NaN's.
x | Coordinates of data points (dimension nr). |
y | Coordinates of data points (dimension nr). |
nr | Number of data points. |
slope | Estimated slope. |
ic | Estimated intercept. |
ssd | Sum of squared distances / nr. |
Definition at line 453 of file llsqwt.c.
|
extern |
Algorithm for weighting the problem that is given to a LLSQ algorithm.
Square roots of weights are used because in LLSQ algorithms the difference w*A-w*b is squared.
N | Matrix A dimension N (nr of parameters). |
M | Matrix A dimension M (nr of samples). |
A | Pointer to matrix A[N][M]; enter NULL to use following matrix format instead. |
a | Pointer to matrix A[N*M]; enter NULL to use previous matrix format instead. |
b | Vector B of length M. |
weight | Weights for each sample (array of length M). |
Definition at line 419 of file bvls.c.
|
extern |
Algorithm for weighting the problem that is given to a LLSQ algorithm.
Square roots of weights are used because in LLSQ algorithms the difference w*A-w*b is squared.
Here user must give squared weights; this makes calculation faster, when this function needs to be called many times.
N | Matrix A dimension N (nr of parameters). |
M | Matrix A dimension M (nr of samples). |
A | Pointer to matrix A[N][M]; enter NULL to use following matrix format instead. |
a | Pointer to matrix A[N*M]; enter NULL to use previous matrix format instead. |
b | Vector B of length M. |
sweight | Squared weights for each sample (array of length M). |
Definition at line 474 of file bvls.c.
|
extern |
Iterative method for linear least-squares fit with errors in both coordinates.
This function is fully based on article [3].
For n data-point pairs (x[i], y[i]) each point has its own weighting factors in (wx[i], wy[i]). This routine finds the values of the parameters m (slope) and c (intercept, ic) that yield the "best-fit" of the model equation Y = mX + c to the data, where X and Y are the predicted or calculated values of the data points.
Weighting factors wx and wy must be assigned as the inverses of the variances or squares of the measurement uncertainties (SDs), i.e. w[i]=1/(sd[i])^2
If true weights are unknown but yet the relative weights are correct, the slope, intercept and residuals (WSS) will be correct. The applied term S/(N-2) makes also the estimate of sigma (sd) of slope less dependent on the scaling of weights. The sigmas are not exact, since only the lowest-order terms in Taylor-series expansion are incorporated; anyhow sigmas are more accurate than the ones based on York algorithm.
One or more data points can be excluded from the fit by setting either x or y weight to 0.
References:
x | coordinates of data points (of dimension n). |
y | coordinates of data points (of dimension n). |
n | number of data points. |
wx | weighting factors in x. |
wy | weighting factors in y. |
tol | allowed tolerance in slope estimation. |
w | work vector (of dimension n); effective weights w[i] are returned in it. |
ic | Estimated intercept. |
slope | Estimated slope. |
nwss | sqrt(WSS)/wsum of the residuals. |
sic | expected sd of intercept at calculated points; If NULL, then not calculated and variable is left unchanged. |
sslope | Expected sd of slope at calculated points; If NULL, then not calculated and variable is left unchanged. |
cx | Estimated data points (X,y); If NULL, then not calculated and variable is left unchanged. |
cy | Estimated data points (x,Y); If NULL, then not calculated and variable is left unchanged. |
Definition at line 48 of file llsqwt.c.
Referenced by best_llsqwt(), and best_logan_reed().
|
extern |
Calculates Logan plot x,y values from the measured input and ROI concentration TACs.
Plot will not include data where: 1) any of the values is not available (NaN), 2) integral is negative at this or later point, 3) divider is too close to zero.
data_nr | Nr of samples. |
i | Array of input concentrations. |
ii | Array of integrals (from zero to sample time) of input concentrations; if reference region input, then remember to consider the frame length. |
c | Array of ROI concentrations. |
ci | Array of ROI integrals (from zero to frame middle time). |
k2 | Reference region k2; set to <=0 if not needed. |
x | Pointer to preallocated memory (at least size dnr) where MTGA plot x values will be written. |
y | Pointer to preallocated memory (at least size dnr) where MTGA plot y values will be written. |
Definition at line 81 of file mtga.c.
Referenced by img_logan().
|
extern |
Calculates the mean and SD of data. Data (y data) may contain NaN's.
x | Data x values |
y | Data y values |
nr | Number of data sample values |
xmean | Calculated x mean |
xsd | Calculated SD of x mean |
ymean | Calculated y mean |
ysd | Calculated SD of y mean |
Definition at line 341 of file pearson.c.
Referenced by dftMeanTAC(), doubleMatchRel(), imgsegmClusterExpand(), least_trimmed_square(), mertwiRandomExponential(), noiseSD4SimulationFromDFT(), and resMean().
|
extern |
Median-based distribution-free estimation of slope and intercept. This method has no need for weighting and is insensitive to outliers. Note that this is not LMS !
Reference (containing reference to the original idea):
x | Coordinates of data points (dimension nr). |
y | Coordinates of data points (dimension nr). |
nr | Number of data points. |
slope | Estimated slope. |
ic | Estimated intercept. |
Definition at line 533 of file llsqwt.c.
|
extern |
|
extern |
Initialize the state vector mt[] inside data struct for Mersenne Twister MT19937 pseudorandom number generator using given array.
Call either this or mertwiInitWithSeed64 before generating random numbers with MT19937 using functions.
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
init_key | The array for initializing keys |
key_length | Length of initialization array init_key[] |
Definition at line 111 of file mertwi.c.
|
extern |
Initialize the state vector mt[] inside data struct for Mersenne Twister MT19937 pseudorandom number generator using given seed.
Call either this or mertwiInitByArray64 before generating random numbers with MT19937 using functions.
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
seed | Seed, for example from mertwiSeed64() |
Definition at line 91 of file mertwi.c.
Referenced by mertwiInitByArray64(), and mertwiRandomInt64().
|
extern |
Generate a 64-bit double precision floating point pseudorandom number in the range of [0,1] with uniform distribution using Mersenne Twister MT19937.
With uniform distribution, the SD=(up-low)/sqrt(12), and CV=(up-low)/(sqrt(3)*(low+up)).
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
Definition at line 218 of file mertwi.c.
Referenced by mertwiRandomBetween(), mertwiRandomExponential(), and mertwiRandomGaussian().
|
extern |
Generate a 64-bit double precision floating point pseudorandom number in the range of [0,1) with uniform distribution using Mersenne Twister MT19937.
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
|
extern |
Generate a 64-bit double precision floating point pseudorandom number in the range of (0,1) with uniform distribution using Mersenne Twister MT19937.
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
|
extern |
Generate a random number on [0, 2^63-1]-interval using Mersenne Twister MT19937.
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
|
extern |
Generate a random number on [0, 2^64-1]-interval using Mersenne Twister MT19937.
mt | Data struct for Mersenne Twister MT19937 pseudorandom number generator |
Definition at line 152 of file mertwi.c.
Referenced by mertwiRandomDouble1(), mertwiRandomDouble2(), mertwiRandomDouble3(), and mertwiRandomInt63().
|
extern |
Make uint32_t seed for pseudorandom number generators.
Uses microseconds from the computer clock and process ID to reduce the chance of getting the same seed for simultaneously executing program threads and instances.
Definition at line 43 of file mertwi.c.
Referenced by mertwiSeed64().
|
extern |
Make uint64_t seed for pseudorandom number generators.
Uses microseconds from the computer clock and process ID to reduce the chance of getting the same seed for simultaneously executing program threads and instances.
Definition at line 72 of file mertwi.c.
|
extern |
Fit a constant (horizontal straight line) to the data with M-estimator.
The algorithm was described in the lecture notes of Nonlinear signal processing course of Tampere university of technology www.cs.tut.fi/~eeroh/nonlin.html
data | Data array |
nr | Number of data values |
iterNr | Number of iterations |
cutoff | cutoff point |
Definition at line 18 of file mestim.c.
|
extern |
Calculates K1/k2 ratio for [O-15]O2 in muscle, based on OER.
OER | Oxygen extraction ratio; assumptions hold only in the range OER=0.2-0.9 |
SaO2 | Arterial oxygen saturation fraction; you can enter DEFAULT_SAO2 |
p50Hb | Half-saturation pressure p50 (kPa) for hemoglobin; you can enter DEFAULT_P50HB |
p50Mb | Half-saturation pressure p50 (kPa) p50 for myoglobin; you can enter DEFAULT_P50MB |
nHb | Hill coefficient n for hemoglobin; you can enter DEFAULT_NHB |
cHb | Hemoglobin concentration in blood (mg/g); you can enter DEFAULT_CHB |
cMb | Myoglobin concentration in muscle (mg/g); you can enter DEFAULT_CMB |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr |
Definition at line 29 of file o2.c.
|
extern |
Calculates the partial pressure of oxygen in muscle, based on OER and K1/k2.
OER | Oxygen extraction ratio; assumptions hold only in the range OER=0.2-0.9 |
K1k2 | K1/k2 ratio for [O-15]O2 |
SaO2 | Arterial oxygen saturation fraction; you can enter DEFAULT_SAO2 |
p50Mb | Half-saturation pressure p50 (kPa) p50 for myoglobin; you can enter DEFAULT_P50MB |
cHb | Hemoglobin concentration in blood (mg/g); you can enter DEFAULT_CHB |
cMb | Myoglobin concentration in muscle (mg/g); you can enter DEFAULT_CMB |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr |
Definition at line 83 of file o2.c.
|
extern |
Check if model parameters have collided with given limits.
If parameter is fixed (equal lower and upper limit) then it is not counted as collision.
par_nr | Nr of parameters |
lower_p | Lower limits |
upper_p | Upper limits |
test_p | Parameters to test |
Definition at line 59 of file constraints.c.
|
extern |
Check that model parameters are within given limits. If not, then compute a penalty factor.
par_nr | Nr of parameters |
lower_p | Lower limits |
upper_p | Upper limits |
test_p | Parameters to test |
accept_p | Pointer to corrected parameters (NULL if not needed) |
penalty | Pointer to variable in which the possible penalty factor will be written; 1 if no penalty, or >1. Set to NULL if not needed. |
Definition at line 15 of file constraints.c.
|
extern |
Calculates the maximum run length between given two arrays of data.
Definition at line 103 of file runs_test.c.
|
extern |
Finds the best regression line to (x,y)-data, leaving points out from the beginning, because Gjedde-Patlak and Logan plots reach linearity at some later phase.
This function applies llsqperp() which is a non-iterative perpendicular line fitting routine.
x | Plot x axis values. |
y | Plot y axis values. |
nr | Nr of plot data points. |
slope | Slope is returned in here. |
ic | Y axis intercept is returned in here. |
ssd | Sum of squared distances / fnr, or NULL if not needed. |
fnr | Number of points in the best fit, or NULL if not needed. |
Definition at line 141 of file mtga.c.
Referenced by img_logan(), and img_patlak().
|
extern |
Calculates the area under the Gaussian probability density function integrated from minus infinity to given value a. For more information search for the Cephes C codes.
a | variable a |
Definition at line 17 of file normaldistr.c.
Referenced by normal_pvalue_1(), normal_pvalue_2(), and runs_test().
|
extern |
Local one-dimensional minimization by bracketing.
_fun | Pointer to the function to be minimized. It must be defined in main program as: double func(double p, void *data); where p is the parameter, and data points to user-defined data structure needed by the function. |
_fundata | Pointer to data that will be passed on to the _fun(). |
x | Parameter initial value. |
xl | Parameter lower limit. |
xu | Parameter upper limit. |
delta | Initial parameter delta. |
tol | Required tolerance. |
maxeval | Maximum number of function evaluations. |
nx | Pointer for optimized parameter value. |
nf | Pointer for function value at optimum; NULL if not needed. |
verbose | Verbose level |
Definition at line 14 of file nlopt1d.c.
|
extern |
Algorithm NNLS (Non-negative least-squares).
Given an m by n matrix A, and an m-vector B, computes an n-vector X, that solves the least squares problem A * X = B , subject to X>=0
Instead of pointers for working space, NULL can be given to let this function to allocate and free the required memory.
a | On entry, a[ 0... N ][ 0 ... M ] contains the M by N matrix A. On exit, a[][] contains the product matrix Q*A, where Q is an m by n orthogonal matrix generated implicitly by this function. |
m | Matrix dimension m |
n | Matrix dimension n |
b | On entry, b[] must contain the m-vector B. On exit, b[] contains Q*B |
x | On exit, x[] will contain the solution vector |
rnorm | On exit, rnorm contains the squared Euclidean norm of the residual vector, R^2. If NULL is given, no rnorm is calculated |
wp | An n-array of working space, wp[]. On exit, wp[] will contain the dual solution vector. wp[i]=0.0 for all i in set p and wp[i]<=0.0 for all i in set z. Can be NULL, which causes this algorithm to allocate memory for it. |
zzp | An m-array of working space, zz[]. Can be NULL, which causes this algorithm to allocate memory for it. |
indexp | An n-array of working space, index[]. Can be NULL, which causes this algorithm to allocate memory for it. |
Definition at line 38 of file nnls.c.
Referenced by fitExpDecayNNLS(), and img_k1_using_ki().
|
extern |
Algorithm for weighting the problem that is given to nnls-algorithm.
Square roots of weights are used because in nnls the difference w*A-w*b is squared.
N | NNLS dimension N (nr of parameters). |
M | NNLS dimension M (nr of samples). |
A | NNLS matrix A. |
b | NNLS vector B. |
weight | Weights for each sample (array of length M). |
Definition at line 254 of file nnls.c.
|
extern |
Algorithm for weighting the problem that is given to nnls-algorithm.
Square roots of weights are used because in nnls the difference w*A-w*b is squared. Here user must give squared weights; this makes calculation faster, when this function needs to be called many times.
N | NNLS dimension N (nr of parameters). |
M | NNLS dimension M (nr of samples). |
A | NNLS matrix A. |
b | NNLS vector B. |
sweight | Squared weights for each sample (array of length M). |
Definition at line 303 of file nnls.c.
|
extern |
Calculates the one-sided p-value for x in relation to the standard normal distribution (that is, the probability that a random variable distributed as N(0, 1) is greater than x).
x | double-precision value x |
Definition at line 60 of file normaldistr.c.
|
extern |
Calculates the two-sided p-value for x in relation to the standard normal distribution.
x | double-precision value x |
Definition at line 43 of file normaldistr.c.
|
extern |
Calculate the number of free parameters.
Model parameters can be fixed by setting lower and upper limit to equal values. This function simply checks the limits for each parameter.
n | Nr of parameters |
pLower | Lower limits (array of length n) |
pUpper | Upper limits (array of length n) |
Definition at line 50 of file aic.c.
|
extern |
Calculates Gjedde-Patlak plot x,y values from the measured input and ROI concentration TACs.
Plot will not include data where: 1) any of the values is not available (NaN), 2) integral is negative at this or later point, 3) divider is too close to zero, 4) plot x value is negative.
data_nr | Nr of samples. |
i | Array of input concentrations. |
ii | Array of integrals (from zero to sample time) of input concentrations; if reference region input, then remember to consider the frame length. |
c | Array of ROI concentrations. |
x | Pointer to preallocated memory (at least size dnr) where MTGA plot x values will be written. |
y | Pointer to preallocated memory (at least size dnr) where MTGA plot y values will be written. |
Definition at line 22 of file mtga.c.
Referenced by img_patlak().
|
extern |
Calculate slope and intercept of a line and Pearson's correlation coefficient.
x | data x values |
y | data y values |
nr | number of data sample values |
k | slope |
kSD | S.D. of slope |
b | y axis intercept |
bSD | S.D. of y axis intercept |
r | Pearson's correlation coefficient, r |
ySD | Residual variance of y values |
Definition at line 14 of file pearson.c.
Referenced by best_logan_regr(), best_pearson(), dft_end_line(), extrapolate_monoexp(), pearson2(), pearson3(), and pearson4().
|
extern |
Calculate slope and intercept of a line and Pearson's correlation coefficient.
Array char is[] specifies whether single (x,y) points are used in the fit.
x | data x values |
y | data y values |
is | Switch values: 0=do not use this point |
nr | number of data sample values |
k | slope |
kSD | S.D. of slope |
b | y axis intercept |
bSD | S.D. of y axis intercept |
r | Pearson's correlation coefficient, r |
ySD | Residual variance of y values |
Definition at line 109 of file pearson.c.
|
extern |
Calculate slope and intercept of a line and Pearson's correlation coefficient.
Data points may contain NaN's.
x | data x values |
y | data y values |
nr | number of data sample values |
k | slope |
kSD | S.D. of slope |
b | y axis intercept |
bSD | S.D. of y axis intercept |
r | Pearson's correlation coefficient, r |
ySD | Residual variance of y values |
Definition at line 154 of file pearson.c.
|
extern |
Calculate slope and intercept of a line and Pearson's correlation coefficient.
Data points may contain NaN's. Fit start and end times are specified.
x | data x values |
y | data y values |
nr | number of data sample values |
start | fit start time |
end | fit end time |
k | slope |
kSD | S.D. of slope |
b | y axis intercept |
bSD | S.D. of y axis intercept |
r | Pearson's correlation coefficient, r |
ySD | Residual variance of y values |
Definition at line 198 of file pearson.c.
|
extern |
Integrate PET TAC data to frame mid times.
Any of output arrays may be set to NULL if that is not needed. Frames must be in ascending time order. Gaps and small overlap are allowed. If x1[0]>0 and x1[0]<=x2[0]-x1[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x1 | frame start times |
x2 | frame end times |
y | avg value during frame |
nr | number of frames |
ie | integrals at frame mid time |
iie | 2nd integrals at frame mid time |
Definition at line 771 of file integr.c.
Referenced by dftInterpolate(), dftInterpolateForIMG(), dftInterpolateInto(), dftReadinput(), img_k1_using_ki(), img_logan(), and img_patlak().
|
extern |
Calculates integrals of PET data at frame end times.
Data does not have to be continuous, but it must be increasing in time. For faster performance in repetitive calls, allocate memory for integral[] even if it is not needed. If x1[0] is >0 AND x1[0] is <=(x2[0]-x1[0]), then the beginning is interpolated from (0, 0) to (x1[0], y1[0]*x1[0]/x) where x is midtime of the first frame.
x1 | Array of frame start times |
x2 | Array of frame end times |
y | Array of y values (avg during each frame) |
nr | Nr of frames |
newyi | Output: integral values at frame end times, or NULL |
newyii | Output: 2nd integral values at frame end times, or NULL |
Definition at line 334 of file integr.c.
Referenced by dftReadReference().
|
extern |
Integrate PET TAC data to frame end times.
Any of output arrays may be set to NULL if that is not needed. Frames must be in ascending time order. Gaps and small overlap are allowed. If x1[0]>0 and x1[0]<=x2[0]-x1[0], then an imaginary line is drawn from (0,0) to (x[0],y[0]).
x1 | frame start times |
x2 | frame end times |
y | avg value during frame |
nr | number of frames |
e | values at frame end time |
ie | integrals at frame end time |
iie | 2nd integrals at frame end time |
Definition at line 905 of file integr.c.
|
extern |
Powell function minimization routine.
p | Initial guess and final set of parameters |
delta | Initial changes for parameters, ==0 if fixed |
parNr | Nr of parameters |
ftol | Fractional tolerance (for WSS); 0<ftol<1 |
iterNr | Max nr of iterations, and nr of required iters |
fret | Function return value (WSS) at minimum |
_fun | Function to minimize (must return the WSS) |
fundata | Pointer to data which is passed on to the function; NULL if not needed |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr |
Definition at line 43 of file powell.c.
Referenced by bootstrap(), and tgo().
|
extern |
Algorithm QR.
Solves a matrix form least square problem min||A x - b|| => A x =~ b (A is m*n matrix, m>=n) using the QR decomposition for overdetermined systems. Based on GNU Scientific Library, edited by Kaisa Liukko.
Instead of pointers for working space, NULL can be given to let this function to allocate and free the required memory.
A | On entry, a[m][n] contains the m by n matrix A. On exit, a[][] contains the QR factorization. |
m | Dimensions of matrix A are a[m][n]. |
n | Dimensions of matrix A are a[m][n]. |
B | B[] is an m-size vector containing the right-hand side vector b. |
X | On exit, x[] will contain the solution vector x (size of n). |
rnorm | On exit, rnorm (pointer to double) contains the squared Euclidean norm of the residual vector (R^2); enter NULL if not needed. |
tau | On exit, tau[] will contain the householder coefficients (size of n); enter NULL, if not needed. |
res | An m-size array of working space, res[]. On output contains residual b - Ax. Enter NULL to let qr() to handle it. |
wws | m*n array of working space. Enter NULL to let qr() to handle it. |
ws | 2m-array of working space. Enter NULL to let qr() to handle it. |
Definition at line 346 of file qr.c.
|
extern |
Factorise a general M x N matrix A into A = Q R , where Q is orthogonal (M x M) and R is upper triangular (M x N).
Q is stored as a packed set of Householder vectors in the strict lower triangular part of the input matrix A and a set of coefficients in vector tau.
R is stored in the diagonal and upper triangle of the input matrix.
The full matrix for Q can be obtained as the product Q = Q_1 Q_2 .. Q_k and it's transform as the product Q^T = Q_k .. Q_2 Q_1 , where k = min(M,N) and Q_i = (I - tau_i * h_i * h_i^T) and where h_i is a Householder vector h_i = [1, A(i+1,i), A(i+2,i), ... , A(M,i)]. This storage scheme is the same as in LAPACK.
NOTICE! The calling program must take care that pointer tau is of size N or M, whichever is smaller.
a | contains coefficient matrix A (m*n) as input and factorisation QR as output. |
M | nr of rows in matrix A. |
N | nr of columns in matrix A. |
tau | Vector for householder coefficients, of length N or M, whichever is smaller. |
cchain | m*n matrix of working space. |
chain | m size array of working space. |
Definition at line 427 of file qr.c.
Referenced by qr().
|
extern |
Find the least squares solution to the overdetermined system
A x = b
for m >= n using the QR factorisation A = Q R. qr_decomp() must be used prior to this function in order to form the QR factorisation of A. Solution is formed in the following order: QR x = b => R x = Q^T b => x = R^-1 (Q^T b)
NOTICE! The calling program must take care that pointers b, x and residual are of the right size.
QR | m*n matrix containing householder vectors of A. |
M | nr of rows in matrix A. |
N | nr of columns in matrix A. |
tau | vector containing householder coefficients tau; length is M or N, whichever is smaller. |
b | Contains m-size vector b of A x = b. |
x | solution vector x of length n. |
residual | residual vector of length m. |
resNorm | norm^2 of the residual vector; enter NULL if not needed. |
cchain | m*n matrix of the working space. |
chain | 2m length array of the working space. |
Definition at line 492 of file qr.c.
Referenced by qr().
|
extern |
Algorithm for weighting the problem that is given to QR algorithm.
Square roots of weights are used because in QR the difference w*A-w*b is squared.
N | Dimensions of matrix A are a[m][n]. |
M | Dimensions of matrix A are a[m][n]; size of vector B is m. |
A | Matrix a[m][n] for QR, contents will be weighted here. |
b | B[] is an m-size vector for QR, contents will be weighted here. |
weight | Pointer to array of size m, which contains sample weights, used here to weight matrix A and vector b. |
ws | m-sized vector for working space; enter NULL to allocate locally. |
Definition at line 576 of file qr.c.
|
extern |
Solve over-determined least-squares problem A x ~ b using successive Householder rotations.
This routine is based on the text and Fortran code in C.L. Lawson and R.J. Hanson, Solving Least Squares Problems, Prentice-Hall, Englewood Cliffs, New Jersey, 1974, and Fortran code by R.L. Parker and P.B. Stark.
m | Number of samples in matrix A and the length of vector b. |
n | Number of parameters in matrix A and the length of vector x. The n must be smaller or equal to m. |
a | Pointer to matrix A; matrix must be given as an n*m array, containing n consecutive m-length vectors. Contents of A are modified in this routine. |
b | Pointer to vector b of length n. Contents of b are modified in this routine. |
x | Pointer to the result vector x of length n. |
r2 | Pointer to a double value, in where the sum of squared residuals is written. |
Definition at line 633 of file qr.c.
Referenced by bvls().
|
extern |
QR least-squares solving routine.
Solves a matrix form least square problem min||A x - b|| => A x =~ b (A is m*n matrix, m>=n) using the QR decomposition for overdetermined systems
mat | Pointer to the row-major matrix (2D array), A[rows][cols]; not modified. |
rhs | Vector b[rows]; modified to contain the computed (fitted) b[], that is, matrix-vector product A*x. |
sol | Solution vector x[cols]; solution x corresponds to the solution of both the modified and original systems A and B. |
rows | Number of rows (samples) in matrix A and length of vector B. |
cols | Number of cols (parameters) in matrix A and length of vector X. |
r2 | Pointer to value where R^2, the difference between original and computed right hand side (b); enter NULL if not needed. |
Definition at line 54 of file qr.c.
|
extern |
Finds the real roots of a*x^2 + b*x + c = 0
a | Input A |
b | Input B |
c | Input C |
m1 | Output: Root 1 |
m2 | Output: Root 2 |
Definition at line 490 of file llsqwt.c.
Referenced by llsqperp().
|
extern |
Fills the double array with random numbers between specified limits. Set seed for random number generator before calling this routine, for example with srand(time(NULL));
nr | Nr of values in double array |
d | Pointer to allocated double array |
low | Lower limit for random values |
up | Upper limit for random values |
type | Distribution: 0=even, 1=square-root transformation |
Definition at line 159 of file gaussdev.c.
|
extern |
Random shuffle: arrange the n elements of array in random order. Only effective if N is much smaller than RAND_MAX.
Definition at line 13 of file shuffle.c.
Referenced by randperm().
|
extern |
Random permutation: given (allocated) array of length n is filled with random permutation of numbers in the range [a:n+a-1]; that is, each number once in random order.
Definition at line 29 of file shuffle.c.
|
extern |
Calculates regression line slope (m) and y axis intercept.
Data (x and y data) may contain NaN's.
x | An array of x axis values. |
y | An array of y axis values. |
n | The number of values in x and y arrays. |
m | Pointer where calculated slope is written. |
c | Pointer where calculated y axis intercept is written. |
Definition at line 387 of file pearson.c.
Referenced by highest_slope(), and highest_slope_after().
|
extern |
Residuals function calculates the nr of runs (a sequence of consecutive positive or negative residuals) and nr of negative and positive residuals for two data arrays.
d1 | double array with data points |
d2 | double array with data points |
Nr | number of data points in both arrays |
Rnr | nr of runs |
Nminus | nr of negative residuals |
Nplus | nr of positive residuals |
Definition at line 56 of file runs_test.c.
Referenced by runs_test().
|
extern |
Runs_test function tests if residuals of two data arrays are independent.
data1 | double array with data points |
data2 | double array with data points |
N | number of data points in both arrays |
alpha | significance level percentage (set to negative value to use default 5%) |
p | calculated probability p; enter NULL if not needed |
Definition at line 12 of file runs_test.c.
|
extern |
Simulates tissue TAC using 1 tissue compartment model plasma TAC, at plasma TAC times.
Memory for ct must be allocated in the calling program.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
ca | Array of arterial activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
ct | Pointer for TAC array to be simulated; must be allocated |
Definition at line 1317 of file simulate.c.
Referenced by bf_srtm(), bfIrr2TCM(), bfRadiowater(), and simDispersion().
|
extern |
Simulates tissue TAC using 2 tissue compartment model (in series) and plasma TAC, at plasma TAC times. In contrary to the common model, kLoss represents a direct loss rate from the 2nd tissue compartment to venous plasma.
Memory for ct must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta and ctb can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
ca | Array of arterial activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
kLoss | Rate constant of the model |
ct | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
Definition at line 500 of file simulate.c.
|
extern |
Simulates tissue TAC using 2 tissue compartment model and plasma TAC, at plasma TAC times, considering also arterial and venous vasculature. The efflux from 2nd tissue compartment (at rate kL) goes directly to blood.
Memory for cpet must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta, ctb, ctab and/or ctvb can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab.,
t | Array of time values |
ca | Array of arterial plasma activities |
cb | Array of arterial blood activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
kL | Rate constant of the model |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab. |
vb | Vascular volume fraction |
fa | Arterial fraction of vascular volume |
cpet | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
ctab | Pointer for arterial TAC in tissue, or NULL |
ctvb | Pointer for venous TAC in tissue, or NULL |
Definition at line 592 of file simulate.c.
|
extern |
Simulates tissue TAC using dual-input tissue compartment model (1-3 compartments in series for tracer1, and 1 compartment for tracer2) at plasma TAC times, considering also contribution of arterial and venous vasculature, but no exchange between compartments for tracer1 and tracer2. The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec. If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab.
t | Array of time values |
ca1 | Array of arterial plasma activities of tracer1 |
ca2 | Array of arterial plasma activities of tracer2 |
cb | Array of arterial blood activities |
nr | Number of values in TACs |
k1 | Rate constant of the model for tracer1 (from plasma to C1) |
k2 | Rate constant of the model for tracer1 (from C1 to plasma) |
k3 | Rate constant of the model for tracer1 (from C1 to C2) |
k4 | Rate constant of the model for tracer1 (from C2 to C1) |
k5 | Rate constant of the model for tracer1 (from C2 to C3) |
k6 | Rate constant of the model for tracer1 (from C3 to C2) |
k1b | Rate constant of the model for tracer2 (from plasma to C4) |
k2b | Rate constant of the model for tracer2 (from C4 to plasma) |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab. |
vb | Vascular volume fraction |
fa | Arterial fraction of vascular volume |
scpet | Pointer for TAC array to be simulated; must be allocated |
sct1 | Pointer for 1st tracer1 compartment TAC, or NULL if not needed |
sct2 | Pointer for 2nd tracer1 compartment TAC, or NULL if not needed |
sct3 | Pointer for 3rd tracer1 compartment TAC, or NULL if not needed |
sct1b | Pointer for 1st tracer2 compartment TAC, or NULL if not needed |
sctab | Pointer for arterial TAC in tissue, or NULL if not needed |
sctvb | Pointer for venous TAC in tissue, or NULL if not needed |
Definition at line 1387 of file simulate.c.
|
extern |
Simulates tissue TAC using 1-3 tissue compartment model (2nd and 3rd compartments in parallel) and plasma TAC, at plasma TAC times.
Memory for ct must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta, ctb and/or ctc can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
ca | Array of arterial activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
k4 | Rate constant of the model |
k5 | Rate constant of the model |
k6 | Rate constant of the model |
ct | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
ctc | Pointer for 3rd compartment TAC to be simulated, or NULL |
Definition at line 136 of file simulate.c.
|
extern |
Simulates tissue TAC using 1-3 tissue compartment model (in series) and plasma TAC, at plasma TAC times.
Memory for ct must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta, ctb and/or ctc can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
ca | Array of arterial activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
k4 | Rate constant of the model |
k5 | Rate constant of the model |
k6 | Rate constant of the model |
ct | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
ctc | Pointer for 3rd compartment TAC to be simulated, or NULL |
Definition at line 27 of file simulate.c.
|
extern |
Simulates tissue TAC using 1-3 tissue compartment model (2nd and 3rd compartments in parallel) and plasma TAC, at plasma TAC times, considering also arterial and venous vasculature.
Memory for cpet must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta, ctb, ctc, ctab and/or ctvb can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab.,
t | Array of time values |
ca | Array of arterial plasma activities |
cb | Array of arterial blood activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
k4 | Rate constant of the model |
k5 | Rate constant of the model |
k6 | Rate constant of the model |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab. |
vb | Vascular volume fraction |
fa | Arterial fraction of vascular volume |
cpet | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
ctc | Pointer for 3rd compartment TAC to be simulated, or NULL |
ctab | Pointer for arterial TAC in tissue, or NULL |
ctvb | Pointer for venous TAC in tissue, or NULL |
Definition at line 373 of file simulate.c.
|
extern |
Simulates tissue TAC using 3 tissue compartmental model with two parallel compartments, and plasma TAC, at plasma TAC sample times, considering also arterial and venous vasculature. The efflux from 3rd tissue compartment (C) goes directly to blood at rate kLoss.
Memory for cpet must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta, ctb, ctab and/or ctvb can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab.,
t | Array of sample times |
ca | Array of arterial plasma activities |
cb | Array of arterial blood activities |
nr | Number of sample values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
k4 | Rate constant of the model |
k5 | Rate constant of the model |
k6 | Rate constant of the model |
kLoss | Rate constant of the model |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab. |
vb | Vascular volume fraction |
fa | Arterial fraction of vascular volume |
cpet | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st tissue compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
ctc | Pointer for 3rd compartment TAC to be simulated, or NULL |
ctab | Pointer for arterial TAC in tissue, or NULL |
ctvb | Pointer for venous TAC in tissue, or NULL |
Definition at line 707 of file simulate.c.
|
extern |
Simulates tissue TAC using 1-3 tissue compartment model (in series) and plasma TAC, at plasma TAC times, considering also arterial and venous vasculature.
Memory for cpet must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cta, ctb, ctc, ctab and/or ctvb can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab.,
t | Array of time values |
ca | Array of arterial plasma activities |
cb | Array of arterial blood activities |
nr | Number of values in TACs |
k1 | Rate constant of the model |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
k4 | Rate constant of the model |
k5 | Rate constant of the model |
k6 | Rate constant of the model |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab. |
vb | Vascular volume fraction |
fa | Arterial fraction of vascular volume |
cpet | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
ctc | Pointer for 3rd compartment TAC to be simulated, or NULL |
ctab | Pointer for arterial TAC in tissue, or NULL |
ctvb | Pointer for venous TAC in tissue, or NULL |
Definition at line 243 of file simulate.c.
|
extern |
Simulates tissue TAC using dual-input tissue compartment model (compartments 2 and 3 in parallel for tracer1, and 1 compartment for tracer2) at plasma TAC sample times, considering also contribution of arterial and venous vasculature, and transfer of tracer1 to tracer2. The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec. If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab. Reference: TPCMOD0001 Appendix C. Tested with program p2t_di -parallel.
t | Array of time values |
ca1 | Array of arterial plasma activities of tracer1 (parent tracer) |
ca2 | Array of arterial plasma activities of tracer2 (metabolite) |
cb | Array of (total) arterial blood activities |
nr | Number of values in TACs |
k1 | Rate constant of the model for tracer1 (from plasma to C1) |
k2 | Rate constant of the model for tracer1 (from C1 to plasma) |
k3 | Rate constant of the model for tracer1 (from C1 to C2) |
k4 | Rate constant of the model for tracer1 (from C2 to C1) |
k5 | Rate constant of the model for tracer1 (from C1 to C3) |
k6 | Rate constant of the model for tracer1 (from C3 to C1) |
k7 | Rate constant of the model for tracer1 (from C3 to plasma) |
km | Rate constant of the model (from tracer1 in C1 to tracer2 in C4) |
k1b | Rate constant of the model for tracer2 (from plasma to C4) |
k2b | Rate constant of the model for tracer2 (from C4 to plasma) |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab |
vb | Vascular volume fraction (0<=Vb<1) |
fa | Arterial fraction of vascular volume (0<=fa<=1) |
scpet | Pointer for TAC array to be simulated; must be allocated |
sct1 | Pointer for 1st tracer1 compartment TAC, or NULL if not needed |
sct2 | Pointer for 2nd tracer1 compartment TAC, or NULL if not needed |
sct3 | Pointer for 3rd tracer1 compartment TAC, or NULL if not needed |
sct1b | Pointer for 1st tracer2 compartment TAC, or NULL if not needed |
sctab | Pointer for arterial TAC in tissue, or NULL if not needed |
sctvb | Pointer for venous TAC in tissue, or NULL if not needed |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr |
Definition at line 1533 of file simulate.c.
|
extern |
Simulates tissue TAC using dual-input tissue compartment model (1-3 compartments in series for tracer1, and 1 compartment for tracer2) at plasma TAC times, considering also contribution of arterial and venous vasculature, and transfer of tracer1 to tracer2. The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec. If blood flow is set to 0, function assumes that f>>k1, and Cvb=Cab. Reference: TPCMOD0001 Appendix B. Tested with program p2t_di -series.
t | Array of time values |
ca1 | Array of arterial plasma activities of tracer1 (parent tracer) |
ca2 | Array of arterial plasma activities of tracer2 (metabolite) |
cb | Array of (total) arterial blood activities |
nr | Number of values in TACs |
k1 | Rate constant of the model for tracer1 (from plasma to C1) |
k2 | Rate constant of the model for tracer1 (from C1 to plasma) |
k3 | Rate constant of the model for tracer1 (from C1 to C2) |
k4 | Rate constant of the model for tracer1 (from C2 to C1) |
k5 | Rate constant of the model for tracer1 (from C2 to C3) |
k6 | Rate constant of the model for tracer1 (from C3 to C2) |
k7 | Rate constant of the model for tracer1 (from C3 to plasma) |
km | Rate constant of the model (from tracer1 in C1 to tracer2 in C4) |
k1b | Rate constant of the model for tracer2 (from plasma to C4) |
k2b | Rate constant of the model for tracer2 (from C4 to plasma) |
f | Blood flow; if 0, function assumes that f>>k1, and Cvb=Cab. |
vb | Vascular volume fraction |
fa | Arterial fraction of vascular volume |
scpet | Pointer for TAC array to be simulated; must be allocated |
sct1 | Pointer for 1st tracer1 compartment TAC, or NULL if not needed |
sct2 | Pointer for 2nd tracer1 compartment TAC, or NULL if not needed |
sct3 | Pointer for 3rd tracer1 compartment TAC, or NULL if not needed |
sct1b | Pointer for 1st tracer2 compartment TAC, or NULL if not needed |
sctab | Pointer for arterial TAC in tissue, or NULL if not needed |
sctvb | Pointer for venous TAC in tissue, or NULL if not needed |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr |
Definition at line 1724 of file simulate.c.
|
extern |
Simulate the effect of dispersion on a time-activity curve.
x | Array of sample times. |
y | Array of sample values, which will be replaced here by dispersion added values. |
n | Nr of samples. |
tau1 | First dispersion time constant (zero if no dispersion); in same time unit as sample times. |
tau2 | 2nd dispersion time constant (zero if no dispersion); in same time unit as sample times. |
tmp | Array for temporary data, for at least n samples; enter NULL to let function to allocate and free the temporary space. |
Definition at line 1911 of file simulate.c.
Referenced by fitEvaltac().
|
extern |
Simulation of TACs of parent tracer, and 1-2 of its metabolites in plasma using Huang's compartmental model.
The units of model parameters must be related to the sample time unit; 1/min and min, or 1/sec and sec.
Pointers to memory for output TACs must be specified, or NULL if TAC is not needed.
t | Input: Sample times (preferably with short intervals) |
ctot | Input: Measured total plasma TAC |
nr | Input: Nr of samples |
k01 | Input: Model parameters |
k12 | Input: Model parameters |
k21 | Input: Model parameters |
k03 | Input: Model parameters |
k34 | Input: Model parameters |
k43 | Input: Model parameters |
c0 | Output: unchanged (parent) tracer TAC |
c1 | Output: TAC of the 1st metabolite |
c3 | Output: TAC of the 2nd metabolite |
Definition at line 1054 of file simulate.c.
|
extern |
Simulate myocardial tissue TAC using Iida's compartment model. Memory for ct must be allocated in the calling program. The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
ci | Input activities |
nr | Number of values in TACs |
k1 | Apparent k1 |
k2 | Apparent k2 |
Vfit | Vfit |
ct | Pointer for TAC array to be simulated; must be allocated |
Definition at line 1253 of file simulate.c.
|
extern |
Simulate tissue and venous blood TACs using dual-input compartment model for [O-15]O2 (one tissue compartment for [O-15]O2, and another tissue compartment for its metabolite [O-15]H2O).
The units of rate constants must be related to the time unit of the data; 1/min and min, or 1/sec and sec.
t | Array of sample times |
ca1 | Array of arterial blood activities of tracer1 ([O-15]O2) |
ca2 | Array of arterial blood activities of tracer2 ([O-15]H2O) |
ca1i | Array of AUC 0-t of arterial tracer1 activities; NULL if not available |
ca2i | Array of AUC 0-t of arterial tracer2 activities; NULL if not available |
n | Nr of samples (array lengths) |
k1a | Rate constant of the model for tracer1 (from blood to C1) |
k2a | Rate constant of the model for tracer1 (from C1 to blood) |
km | Rate constant of the model (from tracer1 in C1 to tracer2 in C2) |
k1b | Rate constant of the model for tracer2 (from blood to C2) |
k2b | Rate constant of the model for tracer2 (from C2 to blood) |
vb | Vascular volume fraction [0-1) |
fa | Arterial fraction of vascular volume [0-1] |
scpet | Pointer for TTAC array to be simulated; allocate in the calling program or set to NULL if not needed |
sct1 | Simulated TAC of tracer1 in tissue; allocate in the calling program or set to NULL if not needed |
sct2 | Simulated TAC of tracer2 in tissue; allocate in the calling program or set to NULL if not needed |
sctab | Total arterial contribution to PET TTAC; allocate in the calling program or set to NULL if not needed |
sctvb1 | Venous tracer1 contribution to PET TAC; allocate in the calling program or set to NULL if not needed |
sctvb2 | Venous tarcer1 contribution to PET TAC; allocate in the calling program or set to NULL if not needed |
scvb1 | Venous BTAC of tracer1; allocate in the calling program or set to NULL if not needed |
scvb2 | Venous BTAC of tracer2; allocate in the calling program or set to NULL if not needed |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr |
Definition at line 1978 of file simulate.c.
|
extern |
Downhill simplex function minimization routine.
Note that if any constraints are required for the parameter values they must be set in the function.
_fun | Pointer to the function to be minimized. It must be defined in main program as: double func(double *p); where p is the parameter array. |
parNr | The number of unknown parameters |
par | This double array contains the minimized parameters. Initial values must be set. |
delta | This double array contains the initial changes to parameters. To fix a parameter, set the corresponding delta to 0. |
maxerr | Maximal error allowed (stopping rule #1) |
maxiter | Maximal nr of iterations allowed (stopping rule #2) |
verbose | Verbose level; if zero, then only warnings are printed into stderr |
Definition at line 27 of file simplex.c.
|
extern |
Simulates tissue TAC using reference tissue compartment model (original) and reference region TAC, at reference region TAC times.
Memory for ct must be allocated in the calling program. To retrieve the separate tissue compartment TACs, pointer to allocated memory for cf and/or cb can be given; if compartmental TACs are not required, NULL pointer can be given instead.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
cr | Reference region activities |
nr | Number of values in TACs |
R1 | Ratio K1/K1' |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
k4 | Rate constant of the model |
ct | Pointer for TAC array to be simulated; must be allocated |
cta | Pointer for 1st compartment TAC to be simulated, or NULL |
ctb | Pointer for 2nd compartment TAC to be simulated, or NULL |
Definition at line 835 of file simulate.c.
|
extern |
Simulates tissue TAC using reference tissue compartment model (simplified) and reference region TAC, at reference region TAC times.
Memory for ct must be allocated in the calling program.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
cr | Reference region activities |
nr | Number of values in TACs |
R1 | Ratio K1/K1' |
k2 | Rate constant of the model |
BP | Binding potential |
ct | Pointer for TAC array to be simulated; must be allocated |
Definition at line 919 of file simulate.c.
|
extern |
Simulate parent tracer TAC using plasma metabolite model TPCMOD0009C, https://www.turkupetcentre.net/reports/tpcmod0009_app_c.pdf
t | Sample times |
ctot | Total plasma TAC |
nr | Sample number |
km | km |
k1m | k1m |
k2m | k2m |
k3m | k3m |
k4m | k4m |
ca | Pointer to array where parent tracer TAC will be written; NULL, if not needed. |
cm | Pointer to array where metabolized tracer TAC will be written; NULL, if not needed. |
Definition at line 1165 of file simulate.c.
|
extern |
Simulates tissue TAC using reference tissue compartment model (transport limited in ref region) and reference region TAC, at reference region TAC times.
Memory for ct must be allocated in the calling program.
The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.
t | Array of time values |
cr | Reference region activities |
nr | Number of values in TACs |
R1 | Ratio K1/K1' |
k2 | Rate constant of the model |
k3 | Rate constant of the model |
ct | Pointer for TAC array to be simulated; must be allocated |
Definition at line 986 of file simulate.c.
|
extern |
Topographical minimization algorithm, that searches the global minimum of a function using clusterization. Calls a local minimization algorithm. Based on an algorithm by Aimo Torn and Sami Viitanen.
lowlim | Lower limits for the parameters. |
uplim | Upper limits for the parameters. |
objf | The object function. |
objfData | Data to objective function; NULL if not needed. |
dim | Dimension = nr of parameters. |
neighNr | Nr of neighbours to investigate; enter large nr relative to samNr (below) if only few local minima are expected. |
fmin | Function value at global minimum. |
gmin | Global minimum = parameter estimates. |
samNr | Nr of points to sample in one iteration; enter larger samNr if nr of iterations (below) is small. |
tgoNr | Nr of TGO iterations; enter 0 to use the default; enter 1 to run TGO just once ie TGO instead of iTGO. Large iteration nr is needed if samNr (above) would otherwise require too much memory. |
verbose | Verbose level; if zero, then nothing is printed into stdout or stderr. |
Definition at line 39 of file tgo.c.
|
extern |
Create randomized parameters for TGO.
p | Pointer to list of TGO points. |
parNr | Nr of parameters in the point. |
sNr | Nr of TGO points. |
low | List of lower limits for each parameter. |
up | List of upper limits for each parameter. |
Definition at line 533 of file tgo.c.
Referenced by tgo().
|
extern |
Create randomized parameters for TGO with square-root transformation, that is, parameter distribution is biased towards low absolute value.
p | Pointer to list of TGO points. |
parNr | Nr of parameters in the point. |
sNr | Nr of TGO points. |
low | List of lower limits for each parameter. |
up | List of upper limits for each parameter. |
Definition at line 564 of file tgo.c.
Referenced by tgo().
|
extern |
Seed for random number generator; declared in gaussdev.c
Seed for random number generator
Definition at line 7 of file gaussdev.c.
Referenced by init_gaussdev().
|
extern |
|
extern |