Patterns in static

Apophenia

Functions
Map or apply a function to a vector or matrix

Functions

apop_dataapop_map (apop_data *in, apop_fn_d *fn_d, apop_fn_v *fn_v, apop_fn_r *fn_r, apop_fn_dp *fn_dp, apop_fn_vp *fn_vp, apop_fn_rp *fn_rp, apop_fn_dpi *fn_dpi, apop_fn_vpi *fn_vpi, apop_fn_rpi *fn_rpi, apop_fn_di *fn_di, apop_fn_vi *fn_vi, apop_fn_ri *fn_ri, void *param, int inplace, char part, int all_pages)
 
gsl_vector * apop_matrix_map (const gsl_matrix *m, double(*fn)(gsl_vector *))
 
void apop_matrix_apply (gsl_matrix *m, void(*fn)(gsl_vector *))
 
gsl_vector * apop_vector_map (const gsl_vector *v, double(*fn)(double))
 
void apop_vector_apply (gsl_vector *v, void(*fn)(double *))
 
gsl_matrix * apop_matrix_map_all (const gsl_matrix *in, double(*fn)(double))
 
void apop_matrix_apply_all (gsl_matrix *in, void(*fn)(double *))
 
double apop_vector_map_sum (const gsl_vector *in, double(*fn)(double))
 
double apop_matrix_map_all_sum (const gsl_matrix *in, double(*fn)(double))
 
double apop_matrix_map_sum (const gsl_matrix *in, double(*fn)(gsl_vector *))
 
double apop_map_sum (apop_data *in, apop_fn_d *fn_d, apop_fn_v *fn_v, apop_fn_r *fn_r, apop_fn_dp *fn_dp, apop_fn_vp *fn_vp, apop_fn_rp *fn_rp, apop_fn_dpi *fn_dpi, apop_fn_vpi *fn_vpi, apop_fn_rpi *fn_rpi, apop_fn_di *fn_di, apop_fn_vi *fn_vi, apop_fn_ri *fn_ri, void *param, char part, int all_pages)
 

Detailed Description

These functions will pull each element of a vector or matrix, or each row of a matrix, and apply a function to the given element. See the data->map/apply section of the outline for many examples.

There are two types, which were developed at different times. The apop_map and apop_map_sum functions use variadic function inputs to cover a lot of different types of process depending on the inputs. Other functions with types in their names, like apop_matrix_map and apop_vector_apply, may be easier to use in some cases. With one exception, they use the same guts, so use whichever type is convenient.

Here are a few technical details of usage:

Function Documentation

apop_data* apop_map ( apop_data in,
apop_fn_d *  fn_d,
apop_fn_v *  fn_v,
apop_fn_r *  fn_r,
apop_fn_dp *  fn_dp,
apop_fn_vp *  fn_vp,
apop_fn_rp *  fn_rp,
apop_fn_dpi *  fn_dpi,
apop_fn_vpi *  fn_vpi,
apop_fn_rpi *  fn_rpi,
apop_fn_di *  fn_di,
apop_fn_vi *  fn_vi,
apop_fn_ri *  fn_ri,
void *  param,
int  inplace,
char  part,
int  all_pages 
)

Apply a function to every element of a data set, matrix or vector; or, apply a vector-taking function to every row or column of a matrix.

There are a lot of options: your function could take any combination of a gsl_vector, a double, an apop_data, a parameter set, and the position of the element in the vector or matrix. As such, the function takes twelve function inputs, one for each combination of vector/matrix, params/no params, index/no index. Fortunately, because this function uses the Designated initializers syntax for inputs, you will specify only one.

For example, here is a function that will cut off each element of the input data to between $(-1, +1)$.

1 double cutoff(double in, void *limit_in){
2  double *limit = limit_in;
3  return GSL_MAX(-*limit, GSL_MIN(*limit, in));
4 }
5 
6 double param = 1;
7 apop_map(your_data, .fn_dp=cutoff, .param=&param, .inplace='y');
Parameters
fn_vA function of the form double your_fn(gsl_vector *in)
fn_dA function of the form double your_fn(double in)
fn_rA function of the form double your_fn(apop_data *in)
fn_vpA function of the form double your_fn(gsl_vector *in, void *param)
fn_dpA function of the form double your_fn(double in, void *param)
fn_rpA function of the form double your_fn(apop_data *in, void *param)
fn_vpiA function of the form double your_fn(gsl_vector *in, void *param, int index)
fn_dpiA function of the form double your_fn(double in, void *param, int index)
fn_rpiA function of the form double your_fn(apop_data *in, void *param, int index)
fn_viA function of the form double your_fn(gsl_vector *in, int index)
fn_diA function of the form double your_fn(double in, int index)
fn_riA function of the form double your_fn(apop_data *in, int index)
inThe input data set. If NULL, I'll return NULL immediately.
paramA pointer to the parameters to be passed to those function forms taking a *param.
partWhich part of the apop_data struct should I use?
'v'==Just the vector
'm'==Every element of the matrix, in turn
'a'==Both 'v' and 'm'
'r'==Apply a function gsl_vector $\to$ double to each row of the matrix
'c'==Apply a function gsl_vector $\to$ double to each column of the matrix
Default is 'a', but notice that I'll ignore a NULL vector or matrix, so if your data set has only a vector (for example), that's what I'll use.
all_pagesIf 'y', then I follow the more pointer to subsequent pages, else I handle only the first page of data. [I abuse this for an internal semaphore, by the way, so your input must always be nonnegative and less than 1,000. Of course, 'y' and 'n' fit these rules fine.] Default: 'n'.
inplaceIf 'n' (the default), generate a new apop_data set for output, which will contain the mapped values (and the names from the original set).
If 'y', modify in place. The double $\to$ double versions, 'v', 'm', and 'a', write to exactly the same location as before. The gsl_vector $\to$ double versions, 'r', and 'c', will write to the vector. Be careful: if you are writing in place and there is already a vector there, then the original vector is lost.
If 'v' (as in void), return NULL. (Default = 'n')
  • The function forms with r in them, like fn_ri, are row-by-row. I'll use Apop_r to get each row in turn, and send it to the function. The first implication is that your function should be expecting a apop_data set with exactly one row in it. The second is that part is ignored: it only makes sense to go row-by-row.
  • If you set inplace='y', then you will be modifying your input data set, row by row; if you set inplace='n', then I will return an apop_data set whose vector element is as long as your data set (i.e., as long as the longest of your text, vector, or matrix parts).
  • If you set apop_opts.thread_count to a value greater than one, I will split the data set into as many chunks as you specify, and process them simultaneously. You need to watch out for the usual hang-ups about multithreaded programming, but if your data is iid, and each row's processing is independent of the others, you should have no problems. Bear in mind that generating threads takes some small overhead, so simple cases like adding a few hundred numbers will actually be slower when threading.
Exceptions
out->error='p'missing or mismatched parts error, such as NULL matrix when you sent a function acting on the matrix element.
double apop_map_sum ( apop_data in,
apop_fn_d *  fn_d,
apop_fn_v *  fn_v,
apop_fn_r *  fn_r,
apop_fn_dp *  fn_dp,
apop_fn_vp *  fn_vp,
apop_fn_rp *  fn_rp,
apop_fn_dpi *  fn_dpi,
apop_fn_vpi *  fn_vpi,
apop_fn_rpi *  fn_rpi,
apop_fn_di *  fn_di,
apop_fn_vi *  fn_vi,
apop_fn_ri *  fn_ri,
void *  param,
char  part,
int  all_pages 
)

A function that effectively calls apop_map and returns the sum of the resulting elements. Thus, this function returns a single double. See the apop_map page for details of the inputs, which are the same here, except that inplace doesn't make sense—this function will always just add up the input function outputs.

See also the map/apply page for details.

  • I don't copy the input data to send to your input function. Therefore, if your function modifies its inputs as a side-effect, your data set will be modified as this function runs.
  • The sum of zero elements is zero, so that is what is returned if the input apop_data set is NULL. If apop_opts.verbose >= 2 I print a warning.
void apop_matrix_apply ( gsl_matrix *  m,
void(*)(gsl_vector *)  fn 
)

Apply a function to every row of a matrix. The function that you input takes in a gsl_vector and returns nothing. apop_matrix_apply will produce a vector view of each row, and send each row to your function.

Parameters
mThe matrix
fnA function of the form void fn(gsl_vector* in)
  • If the matrix is NULL, this is a no-op and returns immediately.
  • See the map/apply page for details.
void apop_matrix_apply_all ( gsl_matrix *  in,
void(*)(double *)  fn 
)

Applies a function to every element in a matrix (as opposed to every row)

Parameters
inThe matrix whose elements will be inputs to the function
fnA function with a form like void f(double *in) which will modify the data at the in-pointer in place.
  • If the matrix is NULL, this is a no-op and returns immediately.
  • See the map/apply page for details.
gsl_vector* apop_matrix_map ( const gsl_matrix *  m,
double(*)(gsl_vector *)  fn 
)

Map a function onto every row of a matrix. The function that you input takes in a gsl_vector and returns a double. apop_matrix_map will produce a vector view of each row, and send each row to your function. It will output a gsl_vector holding your function's output for each row.

Parameters
mThe matrix
fnA function of the form double fn(gsl_vector* in)
Returns
A gsl_vector with the corresponding value for each row.
gsl_matrix* apop_matrix_map_all ( const gsl_matrix *  in,
double(*)(double)  fn 
)

Maps a function to every element in a matrix (as opposed to every row)

Parameters
inThe matrix whose elements will be inputs to the function
fnA function with a form like double f(double in).
Returns
a matrix of the same size as the original, with the function applied.
double apop_matrix_map_all_sum ( const gsl_matrix *  in,
double(*)(double)  fn 
)

Like apop_matrix_map_all, but returns the sum of the resulting mapped function. For example, apop_matrix_map_all_sum(v, isnan) returns the number of elements of m that are NaN.

  • If you input a NULL matrix, I return the sum of zero items: zero.
  • See the map/apply page for details.
double apop_matrix_map_sum ( const gsl_matrix *  in,
double(*)(gsl_vector *)  fn 
)

Like apop_matrix_map, but returns the sum of the resulting mapped vector. For example, let log_like be a function that returns the log likelihood of an input vector; then apop_matrix_map_sum(m, log_like) returns the total log likelihood of the rows of m.

  • If you input a NULL matrix, I return the sum of zero items: zero.
  • See the map/apply page for details.
void apop_vector_apply ( gsl_vector *  v,
void(*)(double *)  fn 
)

Apply a function to every row of a matrix. The function that you input takes in a gsl_vector and returns nothing. apop_apply will send a pointer to each element of your vector to your function.

Parameters
vThe input vector
fnA function of the form void fn(double in)
  • If the vector is NULL, this is a no-op and returns immediately.
  • See the map/apply page for details.
gsl_vector* apop_vector_map ( const gsl_vector *  v,
double(*)(double)  fn 
)

Map a function onto every element of a vector. The function that you input takes in a double and returns a double. apop_apply will send each element to your function, and will output a gsl_vector holding your function's output for each row.

Parameters
vThe input vector
fnA function of the form double fn(double in)
Returns
A gsl_vector (allocated by this function) with the corresponding value for each row.
double apop_vector_map_sum ( const gsl_vector *  in,
double(*)(double)  fn 
)

Like apop_vector_map, but returns the sum of the resulting mapped function. For example, apop_vector_map_sum(v, isnan) returns the number of elements of v that are NaN.

  • If you input a NULL vector, I return the sum of zero items: zero.
  • See the map/apply page for details.

Autogenerated by doxygen on Wed Oct 15 2014 (Debian 0.999b+ds3-2).