Patterns in static

Apophenia

apop.h
1 
3 /* Copyright (c) 2005--2014 by Ben Klemens. Licensed under the GPLv2; see COPYING. */
4 
5 /* Here are the headers for all of apophenia's functions, typedefs, static variables and
6 macros. All of these begin with the apop_ (or Apop_ or APOP_) prefix.
7 
8 There used to be a series of sub-headers, but they never provided any serious
9 benefit. Please use your text editor's word-search feature to find any elements you
10 may be looking for. About a third of the file is comments and doxygen documentation,
11 so syntax highlighting that distinguishes code from comments will also help to make
12 this more navigable.*/
13 
19 #pragma once
20 #ifdef __cplusplus
21 extern "C" {
22 #endif
23 
25 #ifndef _GNU_SOURCE
26 #define _GNU_SOURCE //for asprintf
27 #endif
28 
29 #include <assert.h>
30 #include <signal.h> //raise(SIGTRAP)
31 #include <string.h>
32 #include <gsl/gsl_rng.h>
33 #include <gsl/gsl_matrix.h>
34 
35 
37 
38 /* A means of providing more script-like means of sending arguments to a function.
39 
40 These macros are intended as internal. If you are interested in using this mechanism
41 in out-of-Apophenia work, grep docs/documentation.h for optionaldetails to find notes
42 on how these are used (Doxygen doesn't use that page),
43 */
44 #define apop_varad_head(type, name) type variadic_##name(variadic_type_##name varad_in)
45 
46 #define apop_varad_declare(type, name, ...) \
47  typedef struct { \
48  __VA_ARGS__ ; \
49  } variadic_type_##name; \
50  apop_varad_head(type, name);
51 
52 #define apop_varad_var(name, value) name = varad_in.name ? varad_in.name : (value);
53 #define apop_varad_link(name,...) variadic_##name((variadic_type_##name) {__VA_ARGS__})
54  //End of Doxygen ignore.
56 
57 
59 
62 typedef struct{
63  char *title;
64  char * vector;
65  char ** col;
66  char ** row;
67  char ** text;
68  int colct, rowct, textct;
69 } apop_name;
70 
72 typedef struct apop_data{
73  gsl_vector *vector;
74  gsl_matrix *matrix;
75  apop_name *names;
76  char ***text;
77  size_t textsize[2];
78  gsl_vector *weights;
79  struct apop_data *more;
80  char error;
81 } apop_data;
82 
83 /* Settings groups. For internal use only; see apop_settings.c and
84  settings.h for related machinery. */
85 typedef struct {
86  char name[101];
87  unsigned long name_hash;
88  void *setting_group;
89  void *copy;
90  void *free;
92 
94 typedef struct apop_model apop_model;
95 
98 struct apop_model{
99  char name[101];
100  int vsize, msize1, msize2, dsize;
101  apop_data *data;
102  apop_data *parameters;
103  apop_data *info;
104  void (*estimate)(apop_data * data, apop_model *params);
105  long double (*p)(apop_data *d, apop_model *params);
106  long double (*log_likelihood)(apop_data *d, apop_model *params);
107  long double (*cdf)(apop_data *d, apop_model *params);
108  long double (*constraint)(apop_data *data, apop_model *params);
109  int (*draw)(double *out, gsl_rng* r, apop_model *params);
110  void (*prep)(apop_data *data, apop_model *params);
111  apop_settings_type *settings;
112  void *more;
113  size_t more_size;
114  char error;
115 };
116 
118 typedef struct{
119  int verbose;
120  char stop_on_warning;
121  char output_delimiter[100];
123  char input_delimiters[100];
124  char *db_name_column;
125  char *nan_string;
126  char db_engine;
127  char db_user[101];
128  char db_pass[101];
129  FILE *log_file;
132 #define Autoconf_no_atomics 1
133 
134  #if __STDC_VERSION__ > 201100L && !defined(__STDC_NO_ATOMICS__) && Autoconf_no_atomics==0
135  _Atomic(int) rng_seed;
136  #else
137  int rng_seed;
138  #endif
139  float version;
141 
143 
144 apop_name * apop_name_alloc(void);
145 int apop_name_add(apop_name * n, char const *add_me, char type);
146 void apop_name_free(apop_name * free_me);
147 void apop_name_print(apop_name * n);
148 #ifdef APOP_NO_VARIADIC
149  void apop_name_stack(apop_name * n1, apop_name *nadd, char type1, char typeadd) ;
150 #else
151  void apop_name_stack_base(apop_name * n1, apop_name *nadd, char type1, char typeadd) ;
152  apop_varad_declare(void, apop_name_stack, apop_name * n1; apop_name *nadd; char type1; char typeadd);
153 #define apop_name_stack(...) apop_varad_link(apop_name_stack, __VA_ARGS__)
154 #endif
155 
157 int apop_name_find(const apop_name *n, const char *findme, const char type);
158 
159 void apop_data_add_names_base(apop_data *d, const char type, char const ** names);
160 
180 #define apop_data_add_names(dataset, type, ...) apop_data_add_names_base((dataset), (type), (char const*[]) {__VA_ARGS__, NULL})
181 
182 
192 #define apop_data_free(freeme) (apop_data_free_base(freeme) ? 0 : ((freeme)= NULL))
193 
194 char apop_data_free_base(apop_data *freeme);
195 #ifdef APOP_NO_VARIADIC
196  apop_data * apop_data_alloc(const size_t size1, const size_t size2, const int size3) ;
197 #else
198  apop_data * apop_data_alloc_base(const size_t size1, const size_t size2, const int size3) ;
199  apop_varad_declare(apop_data *, apop_data_alloc, const size_t size1; const size_t size2; const int size3);
200 #define apop_data_alloc(...) apop_varad_link(apop_data_alloc, __VA_ARGS__)
201 #endif
202 
203 #ifdef APOP_NO_VARIADIC
204  apop_data * apop_data_calloc(const size_t size1, const size_t size2, const int size3) ;
205 #else
206  apop_data * apop_data_calloc_base(const size_t size1, const size_t size2, const int size3) ;
207  apop_varad_declare(apop_data *, apop_data_calloc, const size_t size1; const size_t size2; const int size3);
208 #define apop_data_calloc(...) apop_varad_link(apop_data_calloc, __VA_ARGS__)
209 #endif
210 
211 #ifdef APOP_NO_VARIADIC
212  apop_data * apop_data_stack(apop_data *m1, apop_data * m2, char posn, char inplace) ;
213 #else
214  apop_data * apop_data_stack_base(apop_data *m1, apop_data * m2, char posn, char inplace) ;
215  apop_varad_declare(apop_data *, apop_data_stack, apop_data *m1; apop_data * m2; char posn; char inplace);
216 #define apop_data_stack(...) apop_varad_link(apop_data_stack, __VA_ARGS__)
217 #endif
218 
219 apop_data ** apop_data_split(apop_data *in, int splitpoint, char r_or_c);
220 apop_data * apop_data_copy(const apop_data *in);
221 void apop_data_rm_columns(apop_data *d, int *drop);
222 void apop_data_memcpy(apop_data *out, const apop_data *in);
223 #ifdef APOP_NO_VARIADIC
224  double * apop_data_ptr(apop_data *data, int row, int col, const char *rowname, const char *colname, const char *page) ;
225 #else
226  double * apop_data_ptr_base(apop_data *data, int row, int col, const char *rowname, const char *colname, const char *page) ;
227  apop_varad_declare(double *, apop_data_ptr, apop_data *data; int row; int col; const char *rowname; const char *colname; const char *page);
228 #define apop_data_ptr(...) apop_varad_link(apop_data_ptr, __VA_ARGS__)
229 #endif
230 
231 #ifdef APOP_NO_VARIADIC
232  double apop_data_get(const apop_data *data, size_t row, int col, const char *rowname, const char *colname, const char *page) ;
233 #else
234  double apop_data_get_base(const apop_data *data, size_t row, int col, const char *rowname, const char *colname, const char *page) ;
235  apop_varad_declare(double, apop_data_get, const apop_data *data; size_t row; int col; const char *rowname; const char *colname; const char *page);
236 #define apop_data_get(...) apop_varad_link(apop_data_get, __VA_ARGS__)
237 #endif
238 
239 #ifdef APOP_NO_VARIADIC
240  int apop_data_set(apop_data *data, size_t row, int col, const double val, const char *rowname, const char * colname, const char *page) ;
241 #else
242  int apop_data_set_base(apop_data *data, size_t row, int col, const double val, const char *rowname, const char * colname, const char *page) ;
243  apop_varad_declare(int, apop_data_set, apop_data *data; size_t row; int col; const double val; const char *rowname; const char * colname; const char *page);
244 #define apop_data_set(...) apop_varad_link(apop_data_set, __VA_ARGS__)
245 #endif
246 
247 void apop_data_add_named_elmt(apop_data *d, char *name, double val);
248 int apop_text_set(apop_data *in, const size_t row, const size_t col, const char *fmt, ...);
249 apop_data * apop_text_alloc(apop_data *in, const size_t row, const size_t col);
250 void apop_text_free(char ***freeme, int rows, int cols);
251 #ifdef APOP_NO_VARIADIC
252  apop_data * apop_data_transpose(apop_data *in, char transpose_text, char inplace) ;
253 #else
254  apop_data * apop_data_transpose_base(apop_data *in, char transpose_text, char inplace) ;
255  apop_varad_declare(apop_data *, apop_data_transpose, apop_data *in; char transpose_text; char inplace);
256 #define apop_data_transpose(...) apop_varad_link(apop_data_transpose, __VA_ARGS__)
257 #endif
258 
259 gsl_matrix * apop_matrix_realloc(gsl_matrix *m, size_t newheight, size_t newwidth);
260 gsl_vector * apop_vector_realloc(gsl_vector *v, size_t newheight);
261 
262 #define apop_data_prune_columns(in, ...) apop_data_prune_columns_base((in), (char *[]) {__VA_ARGS__, NULL})
263 apop_data* apop_data_prune_columns_base(apop_data *d, char **colnames);
264 
265 #ifdef APOP_NO_VARIADIC
266  apop_data * apop_data_get_page(const apop_data * data, const char * title, const char match) ;
267 #else
268  apop_data * apop_data_get_page_base(const apop_data * data, const char * title, const char match) ;
269  apop_varad_declare(apop_data *, apop_data_get_page, const apop_data * data; const char * title; const char match);
270 #define apop_data_get_page(...) apop_varad_link(apop_data_get_page, __VA_ARGS__)
271 #endif
272 
273 apop_data * apop_data_add_page(apop_data * dataset, apop_data *newpage,const char *title);
274 #ifdef APOP_NO_VARIADIC
275  apop_data* apop_data_rm_page(apop_data * data, const char *title, const char free_p) ;
276 #else
277  apop_data* apop_data_rm_page_base(apop_data * data, const char *title, const char free_p) ;
278  apop_varad_declare(apop_data*, apop_data_rm_page, apop_data * data; const char *title; const char free_p);
279 #define apop_data_rm_page(...) apop_varad_link(apop_data_rm_page, __VA_ARGS__)
280 #endif
281 
282 #ifdef APOP_NO_VARIADIC
283  apop_data * apop_data_rm_rows(apop_data *in, int *drop, int (*do_drop)(apop_data* , void*), void* drop_parameter) ;
284 #else
285  apop_data * apop_data_rm_rows_base(apop_data *in, int *drop, int (*do_drop)(apop_data* , void*), void* drop_parameter) ;
286  apop_varad_declare(apop_data *, apop_data_rm_rows, apop_data *in; int *drop; int (*do_drop)(apop_data* , void*); void* drop_parameter);
287 #define apop_data_rm_rows(...) apop_varad_link(apop_data_rm_rows, __VA_ARGS__)
288 #endif
289 
290 
291 //in apop_asst.c:
292 #ifdef APOP_NO_VARIADIC
293  apop_data * apop_model_draws(apop_model *model, int count, apop_data *draws) ;
294 #else
295  apop_data * apop_model_draws_base(apop_model *model, int count, apop_data *draws) ;
296  apop_varad_declare(apop_data *, apop_model_draws, apop_model *model; int count; apop_data *draws);
297 #define apop_model_draws(...) apop_varad_link(apop_model_draws, __VA_ARGS__)
298 #endif
299 
300 
301 
302 /* Convenience functions to convert among vectors (gsl_vector), matrices (gsl_matrix),
303  arrays (double **), and database tables */
304 
305 //From vector
306 gsl_vector *apop_vector_copy(const gsl_vector *in);
307 #ifdef APOP_NO_VARIADIC
308  gsl_matrix * apop_vector_to_matrix(const gsl_vector *in, char row_col) ;
309 #else
310  gsl_matrix * apop_vector_to_matrix_base(const gsl_vector *in, char row_col) ;
311  apop_varad_declare(gsl_matrix *, apop_vector_to_matrix, const gsl_vector *in; char row_col);
312 #define apop_vector_to_matrix(...) apop_varad_link(apop_vector_to_matrix, __VA_ARGS__)
313 #endif
314 
315 
316 //From matrix
317 gsl_matrix *apop_matrix_copy(const gsl_matrix *in);
318 #ifdef APOP_NO_VARIADIC
319  apop_data *apop_db_to_crosstab(char const*tabname, char const*row, char const*col, char const*data, char is_aggregate) ;
320 #else
321  apop_data * apop_db_to_crosstab_base(char const*tabname, char const*row, char const*col, char const*data, char is_aggregate) ;
322  apop_varad_declare(apop_data *, apop_db_to_crosstab, char const*tabname; char const*row; char const*col; char const*data; char is_aggregate);
323 #define apop_db_to_crosstab(...) apop_varad_link(apop_db_to_crosstab, __VA_ARGS__)
324 #endif
325 
326 
327 //From array
328 #ifdef APOP_NO_VARIADIC
329  gsl_vector * apop_array_to_vector(double *in, int size) ;
330 #else
331  gsl_vector * apop_array_to_vector_base(double *in, int size) ;
332  apop_varad_declare(gsl_vector *, apop_array_to_vector, double *in; int size);
333 #define apop_array_to_vector(...) apop_varad_link(apop_array_to_vector, __VA_ARGS__)
334 #endif
335  //Deprecated
337 #define apop_text_add apop_text_set
338 #define apop_line_to_vector apop_array_to_vector
339 
341 //From text
342 #ifdef APOP_NO_VARIADIC
343  apop_data * apop_text_to_data(char const *text_file, int has_row_names, int has_col_names, int const *field_ends, char const *delimiters) ;
344 #else
345  apop_data * apop_text_to_data_base(char const *text_file, int has_row_names, int has_col_names, int const *field_ends, char const *delimiters) ;
346  apop_varad_declare(apop_data *, apop_text_to_data, char const *text_file; int has_row_names; int has_col_names; int const *field_ends; char const *delimiters);
347 #define apop_text_to_data(...) apop_varad_link(apop_text_to_data, __VA_ARGS__)
348 #endif
349 
350 #ifdef APOP_NO_VARIADIC
351  int apop_text_to_db(char const *text_file, char *tabname, int has_row_names, int has_col_names, char **field_names, int const *field_ends, apop_data *field_params, char *table_params, char const *delimiters, char if_table_exists) ;
352 #else
353  int apop_text_to_db_base(char const *text_file, char *tabname, int has_row_names, int has_col_names, char **field_names, int const *field_ends, apop_data *field_params, char *table_params, char const *delimiters, char if_table_exists) ;
354  apop_varad_declare(int, apop_text_to_db, char const *text_file; char *tabname; int has_row_names; int has_col_names; char **field_names; int const *field_ends; apop_data *field_params; char *table_params; char const *delimiters; char if_table_exists);
355 #define apop_text_to_db(...) apop_varad_link(apop_text_to_db, __VA_ARGS__)
356 #endif
357 
358 
359 //rank data
361 #ifdef APOP_NO_VARIADIC
362  apop_data *apop_data_rank_compress (apop_data *in, int min_bins) ;
363 #else
364  apop_data * apop_data_rank_compress_base(apop_data *in, int min_bins) ;
365  apop_varad_declare(apop_data *, apop_data_rank_compress, apop_data *in; int min_bins);
366 #define apop_data_rank_compress(...) apop_varad_link(apop_data_rank_compress, __VA_ARGS__)
367 #endif
368 
369 
370 //From crosstabs
371 void apop_crosstab_to_db(apop_data *in, char *tabname, char *row_col_name,
372  char *col_col_name, char *data_col_name);
373 
374 //packing data into a vector
375 #ifdef APOP_NO_VARIADIC
376  gsl_vector * apop_data_pack(const apop_data *in, gsl_vector *out, char more_pages, char use_info_pages) ;
377 #else
378  gsl_vector * apop_data_pack_base(const apop_data *in, gsl_vector *out, char more_pages, char use_info_pages) ;
379  apop_varad_declare(gsl_vector *, apop_data_pack, const apop_data *in; gsl_vector *out; char more_pages; char use_info_pages);
380 #define apop_data_pack(...) apop_varad_link(apop_data_pack, __VA_ARGS__)
381 #endif
382 
383 #ifdef APOP_NO_VARIADIC
384  void apop_data_unpack(const gsl_vector *in, apop_data *d, char use_info_pages) ;
385 #else
386  void apop_data_unpack_base(const gsl_vector *in, apop_data *d, char use_info_pages) ;
387  apop_varad_declare(void, apop_data_unpack, const gsl_vector *in; apop_data *d; char use_info_pages);
388 #define apop_data_unpack(...) apop_varad_link(apop_data_unpack, __VA_ARGS__)
389 #endif
390 
391 
392 #define apop_vector_fill(avfin, ...) apop_vector_fill_base((avfin), (double []) {__VA_ARGS__})
393 #define apop_data_fill(adfin, ...) apop_data_fill_base((adfin), (double []) {__VA_ARGS__})
394 #define apop_text_fill(dataset, ...) apop_text_fill_base((dataset), (char* []) {__VA_ARGS__, NULL})
395 #define apop_data_falloc(sizes, ...) apop_data_fill(apop_data_alloc sizes, __VA_ARGS__)
396 
397 apop_data *apop_data_fill_base(apop_data *in, double []);
398 gsl_vector *apop_vector_fill_base(gsl_vector *in, double []);
399 apop_data *apop_text_fill_base(apop_data *data, char* text[]);
400 
402 
403 extern apop_model *apop_beta;
404 extern apop_model *apop_bernoulli;
405 extern apop_model *apop_binomial;
406 extern apop_model *apop_chi_squared;
407 extern apop_model *apop_dirichlet;
409 extern apop_model *apop_f_distribution;
410 extern apop_model *apop_gamma;
412 extern apop_model *apop_iv;
414 extern apop_model *apop_loess;
415 extern apop_model *apop_logit;
416 extern apop_model *apop_lognormal;
419 extern apop_model *apop_normal;
420 extern apop_model *apop_ols;
421 extern apop_model *apop_pmf;
422 extern apop_model *apop_poisson;
423 extern apop_model *apop_probit;
425 extern apop_model *apop_uniform;
426 //extern apop_model *apop_wishart;
427 extern apop_model *apop_wls;
428 extern apop_model *apop_yule;
429 extern apop_model *apop_zipf;
430 
431 //model transformations
433 extern apop_model *apop_composition;
435 extern apop_model *apop_mixture;
436 extern apop_model *apop_cross;
437 
439 #define apop_gaussian apop_normal
440 #define apop_OLS apop_ols
441 #define apop_PMF apop_pmf
442 #define apop_F_distribution apop_f_distribution
443 #define apop_IV apop_iv
444 
445 
446 void apop_model_free (apop_model * free_me);
447 #ifdef APOP_NO_VARIADIC
448  void apop_model_print (apop_model * model, FILE *output_pipe) ;
449 #else
450  void apop_model_print_base(apop_model * model, FILE *output_pipe) ;
451  apop_varad_declare(void, apop_model_print, apop_model * model; FILE *output_pipe);
452 #define apop_model_print(...) apop_varad_link(apop_model_print, __VA_ARGS__)
453 #endif
454 
455 void apop_model_show (apop_model * print_me); //deprecated
456 apop_model * apop_model_copy(apop_model *in); //in apop_model.c
458 
460 void apop_score(apop_data *d, gsl_vector *out, apop_model *m);
462 double apop_p(apop_data *d, apop_model *m);
463 double apop_cdf(apop_data *d, apop_model *m);
464 int apop_draw(double *out, gsl_rng *r, apop_model *m);
465 void apop_prep(apop_data *d, apop_model *m);
468 
469 apop_model *apop_beta_from_mean_var(double m, double v); //in apop_beta.c
470 
471 #define apop_model_set_parameters(in, ...) apop_model_set_parameters_base((in), (double []) {__VA_ARGS__})
472 apop_model *apop_model_set_parameters_base(apop_model *in, double ap[]);
473 
474 //apop_mixture.c
480 #define apop_model_mixture(...) apop_model_mixture_base((apop_model *[]){__VA_ARGS__, NULL})
481 apop_model *apop_model_mixture_base(apop_model **inlist);
482 
483 //transform/apop_cross.c.
484 apop_model *apop_model_cross_base(apop_model *mlist[]);
485 #define apop_model_cross(...) apop_model_cross_base((apop_model *[]){__VA_ARGS__, NULL})
486 
488 
489  //The variadic versions, with lots of options to input extra parameters to the
490  //function being mapped/applied
491 #ifdef APOP_NO_VARIADIC
492  apop_data * apop_map(apop_data *in, double (*fn_d)(double), double (*fn_v)(gsl_vector*),
493  double (*fn_r)(apop_data *), double (*fn_dp)(double, void *), double (*fn_vp)(gsl_vector*, void *),
494  double (*fn_rp)(apop_data *, void *), double (*fn_dpi)(double, void *, int),
495  double (*fn_vpi)(gsl_vector*, void *, int), double (*fn_rpi)(apop_data*, void *, int),
496  double (*fn_di)(double, int), double (*fn_vi)(gsl_vector*, int), double (*fn_ri)(apop_data*, int),
497  void *param, int inplace, char part, int all_pages) ;
498 #else
499  apop_data * apop_map_base(apop_data *in, double (*fn_d)(double), double (*fn_v)(gsl_vector*),
500  double (*fn_r)(apop_data *), double (*fn_dp)(double, void *), double (*fn_vp)(gsl_vector*, void *),
501  double (*fn_rp)(apop_data *, void *), double (*fn_dpi)(double, void *, int),
502  double (*fn_vpi)(gsl_vector*, void *, int), double (*fn_rpi)(apop_data*, void *, int),
503  double (*fn_di)(double, int), double (*fn_vi)(gsl_vector*, int), double (*fn_ri)(apop_data*, int),
504  void *param, int inplace, char part, int all_pages) ;
505  apop_varad_declare(apop_data *, apop_map, apop_data *in; double (*fn_d)(double); double (*fn_v)(gsl_vector*);
506  double (*fn_r)(apop_data *); double (*fn_dp)(double, void *); double (*fn_vp)(gsl_vector*, void *);
507  double (*fn_rp)(apop_data *, void *); double (*fn_dpi)(double, void *, int);
508  double (*fn_vpi)(gsl_vector*, void *, int); double (*fn_rpi)(apop_data*, void *, int);
509  double (*fn_di)(double, int); double (*fn_vi)(gsl_vector*, int); double (*fn_ri)(apop_data*, int);
510  void *param; int inplace; char part; int all_pages);
511 #define apop_map(...) apop_varad_link(apop_map, __VA_ARGS__)
512 #endif
513 
514 #ifdef APOP_NO_VARIADIC
515  double apop_map_sum(apop_data *in, double (*fn_d)(double), double (*fn_v)(gsl_vector*),
516  double (*fn_r)(apop_data *), double (*fn_dp)(double, void *), double (*fn_vp)(gsl_vector*, void *),
517  double (*fn_rp)(apop_data *, void *), double (*fn_dpi)(double, void *, int),
518  double (*fn_vpi)(gsl_vector*, void *, int), double (*fn_rpi)(apop_data*, void *, int),
519  double (*fn_di)(double, int), double (*fn_vi)(gsl_vector*, int), double (*fn_ri)(apop_data*, int),
520  void *param, char part, int all_pages) ;
521 #else
522  double apop_map_sum_base(apop_data *in, double (*fn_d)(double), double (*fn_v)(gsl_vector*),
523  double (*fn_r)(apop_data *), double (*fn_dp)(double, void *), double (*fn_vp)(gsl_vector*, void *),
524  double (*fn_rp)(apop_data *, void *), double (*fn_dpi)(double, void *, int),
525  double (*fn_vpi)(gsl_vector*, void *, int), double (*fn_rpi)(apop_data*, void *, int),
526  double (*fn_di)(double, int), double (*fn_vi)(gsl_vector*, int), double (*fn_ri)(apop_data*, int),
527  void *param, char part, int all_pages) ;
528  apop_varad_declare(double, apop_map_sum, apop_data *in; double (*fn_d)(double); double (*fn_v)(gsl_vector*);
529  double (*fn_r)(apop_data *); double (*fn_dp)(double, void *); double (*fn_vp)(gsl_vector*, void *);
530  double (*fn_rp)(apop_data *, void *); double (*fn_dpi)(double, void *, int);
531  double (*fn_vpi)(gsl_vector*, void *, int); double (*fn_rpi)(apop_data*, void *, int);
532  double (*fn_di)(double, int); double (*fn_vi)(gsl_vector*, int); double (*fn_ri)(apop_data*, int);
533  void *param; char part; int all_pages);
534 #define apop_map_sum(...) apop_varad_link(apop_map_sum, __VA_ARGS__)
535 #endif
536 
537 
538  //the specific-to-a-type versions, quicker and easier when appropriate.
539 gsl_vector *apop_matrix_map(const gsl_matrix *m, double (*fn)(gsl_vector*));
540 gsl_vector *apop_vector_map(const gsl_vector *v, double (*fn)(double));
541 void apop_matrix_apply(gsl_matrix *m, void (*fn)(gsl_vector*));
542 void apop_vector_apply(gsl_vector *v, void (*fn)(double*));
543 gsl_matrix * apop_matrix_map_all(const gsl_matrix *in, double (*fn)(double));
544 void apop_matrix_apply_all(gsl_matrix *in, void (*fn)(double *));
545 
546 double apop_vector_map_sum(const gsl_vector *in, double(*fn)(double));
547 double apop_matrix_map_sum(const gsl_matrix *in, double (*fn)(gsl_vector*));
548 double apop_matrix_map_all_sum(const gsl_matrix *in, double (*fn)(double));
549 
550 
551  // Some output routines
552 #ifdef APOP_NO_VARIADIC
553  void apop_matrix_print(const gsl_matrix *data, char const *output_name, FILE *output_pipe, char output_type, char output_append) ;
554 #else
555  void apop_matrix_print_base(const gsl_matrix *data, char const *output_name, FILE *output_pipe, char output_type, char output_append) ;
556  apop_varad_declare(void, apop_matrix_print, const gsl_matrix *data; char const *output_name; FILE *output_pipe; char output_type; char output_append);
557 #define apop_matrix_print(...) apop_varad_link(apop_matrix_print, __VA_ARGS__)
558 #endif
559 
560 #ifdef APOP_NO_VARIADIC
561  void apop_vector_print(gsl_vector *data, char const *output_name, FILE *output_pipe, char output_type, char output_append) ;
562 #else
563  void apop_vector_print_base(gsl_vector *data, char const *output_name, FILE *output_pipe, char output_type, char output_append) ;
564  apop_varad_declare(void, apop_vector_print, gsl_vector *data; char const *output_name; FILE *output_pipe; char output_type; char output_append);
565 #define apop_vector_print(...) apop_varad_link(apop_vector_print, __VA_ARGS__)
566 #endif
567 
568 #ifdef APOP_NO_VARIADIC
569  void apop_data_print(const apop_data *data, char const *output_name, FILE *output_pipe, char output_type, char output_append) ;
570 #else
571  void apop_data_print_base(const apop_data *data, char const *output_name, FILE *output_pipe, char output_type, char output_append) ;
572  apop_varad_declare(void, apop_data_print, const apop_data *data; char const *output_name; FILE *output_pipe; char output_type; char output_append);
573 #define apop_data_print(...) apop_varad_link(apop_data_print, __VA_ARGS__)
574 #endif
575 
576 
577 void apop_matrix_show(const gsl_matrix *data);
578 void apop_vector_show(const gsl_vector *data);
579 void apop_data_show(const apop_data *data);
580 
581 
582  //statistics
583 #ifdef APOP_NO_VARIADIC
584  double apop_vector_mean(gsl_vector const *v, gsl_vector const *weights);
585 #else
586  double apop_vector_mean_base(gsl_vector const *v, gsl_vector const *weights);
587  apop_varad_declare(double, apop_vector_mean, gsl_vector const *v; gsl_vector const *weights);
588 #define apop_vector_mean(...) apop_varad_link(apop_vector_mean, __VA_ARGS__)
589 #endif
590 
591 #ifdef APOP_NO_VARIADIC
592  double apop_vector_var(gsl_vector const *v, gsl_vector const *weights);
593 #else
594  double apop_vector_var_base(gsl_vector const *v, gsl_vector const *weights);
595  apop_varad_declare(double, apop_vector_var, gsl_vector const *v; gsl_vector const *weights);
596 #define apop_vector_var(...) apop_varad_link(apop_vector_var, __VA_ARGS__)
597 #endif
598 
599 #ifdef APOP_NO_VARIADIC
600  double apop_vector_skew_pop(gsl_vector const *v, gsl_vector const *weights);
601 #else
602  double apop_vector_skew_pop_base(gsl_vector const *v, gsl_vector const *weights);
603  apop_varad_declare(double, apop_vector_skew_pop, gsl_vector const *v; gsl_vector const *weights);
604 #define apop_vector_skew_pop(...) apop_varad_link(apop_vector_skew_pop, __VA_ARGS__)
605 #endif
606 
607 #ifdef APOP_NO_VARIADIC
608  double apop_vector_kurtosis_pop(gsl_vector const *v, gsl_vector const *weights);
609 #else
610  double apop_vector_kurtosis_pop_base(gsl_vector const *v, gsl_vector const *weights);
611  apop_varad_declare(double, apop_vector_kurtosis_pop, gsl_vector const *v; gsl_vector const *weights);
612 #define apop_vector_kurtosis_pop(...) apop_varad_link(apop_vector_kurtosis_pop, __VA_ARGS__)
613 #endif
614 
615 #ifdef APOP_NO_VARIADIC
616  double apop_vector_cov(gsl_vector const *v1, gsl_vector const *v2,
617  gsl_vector const *weights);
618 #else
619  double apop_vector_cov_base(gsl_vector const *v1, gsl_vector const *v2,
620  gsl_vector const *weights);
621  apop_varad_declare(double, apop_vector_cov, gsl_vector const *v1; gsl_vector const *v2;
622  gsl_vector const *weights);
623 #define apop_vector_cov(...) apop_varad_link(apop_vector_cov, __VA_ARGS__)
624 #endif
625 
626 
627 #ifdef APOP_NO_VARIADIC
628  double apop_vector_distance(const gsl_vector *ina, const gsl_vector *inb, const char metric, const double norm) ;
629 #else
630  double apop_vector_distance_base(const gsl_vector *ina, const gsl_vector *inb, const char metric, const double norm) ;
631  apop_varad_declare(double, apop_vector_distance, const gsl_vector *ina; const gsl_vector *inb; const char metric; const double norm);
632 #define apop_vector_distance(...) apop_varad_link(apop_vector_distance, __VA_ARGS__)
633 #endif
634 
635 
636 #ifdef APOP_NO_VARIADIC
637  void apop_vector_normalize(gsl_vector *in, gsl_vector **out, const char normalization_type) ;
638 #else
639  void apop_vector_normalize_base(gsl_vector *in, gsl_vector **out, const char normalization_type) ;
640  apop_varad_declare(void, apop_vector_normalize, gsl_vector *in; gsl_vector **out; const char normalization_type);
641 #define apop_vector_normalize(...) apop_varad_link(apop_vector_normalize, __VA_ARGS__)
642 #endif
643 
644 
645 apop_data * apop_data_covariance(const apop_data *in);
646 apop_data * apop_data_correlation(const apop_data *in);
647 long double apop_vector_entropy(gsl_vector *in);
648 long double apop_matrix_sum(const gsl_matrix *m);
649 double apop_matrix_mean(const gsl_matrix *data);
650 void apop_matrix_mean_and_var(const gsl_matrix *data, double *mean, double *var);
651 apop_data * apop_data_summarize(apop_data *data);
652 #ifdef APOP_NO_VARIADIC
653  double * apop_vector_percentiles(gsl_vector *data, char rounding) ;
654 #else
655  double * apop_vector_percentiles_base(gsl_vector *data, char rounding) ;
656  apop_varad_declare(double *, apop_vector_percentiles, gsl_vector *data; char rounding);
657 #define apop_vector_percentiles(...) apop_varad_link(apop_vector_percentiles, __VA_ARGS__)
658 #endif
659 
660 
661 apop_data *apop_test_fisher_exact(apop_data *intab); //in apop_fisher.c
662 
663 //from apop_t_f_chi.c:
664 #ifdef APOP_NO_VARIADIC
665  int apop_matrix_is_positive_semidefinite(gsl_matrix *m, char semi) ;
666 #else
667  int apop_matrix_is_positive_semidefinite_base(gsl_matrix *m, char semi) ;
668  apop_varad_declare(int, apop_matrix_is_positive_semidefinite, gsl_matrix *m; char semi);
669 #define apop_matrix_is_positive_semidefinite(...) apop_varad_link(apop_matrix_is_positive_semidefinite, __VA_ARGS__)
670 #endif
671 
672 double apop_matrix_to_positive_semidefinite(gsl_matrix *m);
673 long double apop_multivariate_gamma(double a, int p);
674 long double apop_multivariate_lngamma(double a, int p);
675 
676 //apop_tests.c
677 apop_data * apop_t_test(gsl_vector *a, gsl_vector *b);
678 apop_data * apop_paired_t_test(gsl_vector *a, gsl_vector *b);
679 #ifdef APOP_NO_VARIADIC
680  apop_data* apop_anova(char *table, char *data, char *grouping1, char *grouping2) ;
681 #else
682  apop_data* apop_anova_base(char *table, char *data, char *grouping1, char *grouping2) ;
683  apop_varad_declare(apop_data*, apop_anova, char *table; char *data; char *grouping1; char *grouping2);
684 #define apop_anova(...) apop_varad_link(apop_anova, __VA_ARGS__)
685 #endif
686 
687 #define apop_ANOVA apop_anova
688 #ifdef APOP_NO_VARIADIC
689  apop_data * apop_f_test (apop_model *est, apop_data *contrast) ;
690 #else
691  apop_data * apop_f_test_base(apop_model *est, apop_data *contrast) ;
692  apop_varad_declare(apop_data *, apop_f_test, apop_model *est; apop_data *contrast);
693 #define apop_f_test(...) apop_varad_link(apop_f_test, __VA_ARGS__)
694 #endif
695 
696 #define apop_F_test apop_f_test
697 
698 //from the regression code:
699 #define apop_estimate_r_squared(in) apop_estimate_coefficient_of_determination(in)
700 
701 apop_data * apop_text_unique_elements(const apop_data *d, size_t col);
702 gsl_vector * apop_vector_unique_elements(const gsl_vector *v);
703 #ifdef APOP_NO_VARIADIC
704  apop_data * apop_data_to_factors(apop_data *data, char intype, int incol, int outcol) ;
705 #else
706  apop_data * apop_data_to_factors_base(apop_data *data, char intype, int incol, int outcol) ;
707  apop_varad_declare(apop_data *, apop_data_to_factors, apop_data *data; char intype; int incol; int outcol);
708 #define apop_data_to_factors(...) apop_varad_link(apop_data_to_factors, __VA_ARGS__)
709 #endif
710 
711 #ifdef APOP_NO_VARIADIC
712  apop_data * apop_data_get_factor_names(apop_data *data, int col, char type) ;
713 #else
714  apop_data * apop_data_get_factor_names_base(apop_data *data, int col, char type) ;
715  apop_varad_declare(apop_data *, apop_data_get_factor_names, apop_data *data; int col; char type);
716 #define apop_data_get_factor_names(...) apop_varad_link(apop_data_get_factor_names, __VA_ARGS__)
717 #endif
718 
719 
720 #ifdef APOP_NO_VARIADIC
721  apop_data * apop_data_to_dummies(apop_data *d, int col, char type, int keep_first, char append, char remove) ;
722 #else
723  apop_data * apop_data_to_dummies_base(apop_data *d, int col, char type, int keep_first, char append, char remove) ;
724  apop_varad_declare(apop_data *, apop_data_to_dummies, apop_data *d; int col; char type; int keep_first; char append; char remove);
725 #define apop_data_to_dummies(...) apop_varad_link(apop_data_to_dummies, __VA_ARGS__)
726 #endif
727 
728 
729 #ifdef APOP_NO_VARIADIC
730  long double apop_model_entropy(apop_model *in, int draws) ;
731 #else
732  long double apop_model_entropy_base(apop_model *in, int draws) ;
733  apop_varad_declare(long double, apop_model_entropy, apop_model *in; int draws);
734 #define apop_model_entropy(...) apop_varad_link(apop_model_entropy, __VA_ARGS__)
735 #endif
736 
737 #ifdef APOP_NO_VARIADIC
738  long double apop_kl_divergence(apop_model *from, apop_model *to, int draw_ct, gsl_rng *rng) ;
739 #else
740  long double apop_kl_divergence_base(apop_model *from, apop_model *to, int draw_ct, gsl_rng *rng) ;
741  apop_varad_declare(long double, apop_kl_divergence, apop_model *from; apop_model *to; int draw_ct; gsl_rng *rng);
742 #define apop_kl_divergence(...) apop_varad_link(apop_kl_divergence, __VA_ARGS__)
743 #endif
744 
745 
747 void apop_estimate_parameter_tests (apop_model *est);
748 
749 //Bootstrapping & RNG
750 apop_data * apop_jackknife_cov(apop_data *data, apop_model *model);
751 #ifdef APOP_NO_VARIADIC
752  apop_data * apop_bootstrap_cov(apop_data *data, apop_model *model, gsl_rng* rng, int iterations, char keep_boots, char ignore_nans, apop_data **boot_store) ;
753 #else
754  apop_data * apop_bootstrap_cov_base(apop_data *data, apop_model *model, gsl_rng* rng, int iterations, char keep_boots, char ignore_nans, apop_data **boot_store) ;
755  apop_varad_declare(apop_data *, apop_bootstrap_cov, apop_data *data; apop_model *model; gsl_rng* rng; int iterations; char keep_boots; char ignore_nans; apop_data **boot_store);
756 #define apop_bootstrap_cov(...) apop_varad_link(apop_bootstrap_cov, __VA_ARGS__)
757 #endif
758 
759 gsl_rng *apop_rng_alloc(int seed);
760 double apop_rng_GHgB3(gsl_rng * r, double* a); //in apop_asst.c
761 
762 #define apop_rng_get_thread(thread_in) apop_rng_get_thread_base(#thread_in[0]=='\0' ? -1: (thread_in+0))
763 gsl_rng *apop_rng_get_thread_base(int thread);
764 
765 int apop_arms_draw (double *out, gsl_rng *r, apop_model *m);
766 
767 
768  // maximum likelihod estimation related functions
769 
770 #ifdef APOP_NO_VARIADIC
771  gsl_vector * apop_numerical_gradient(apop_data * data, apop_model* model, double delta) ;
772 #else
773  gsl_vector * apop_numerical_gradient_base(apop_data * data, apop_model* model, double delta) ;
774  apop_varad_declare(gsl_vector *, apop_numerical_gradient, apop_data * data; apop_model* model; double delta);
775 #define apop_numerical_gradient(...) apop_varad_link(apop_numerical_gradient, __VA_ARGS__)
776 #endif
777 
778 #ifdef APOP_NO_VARIADIC
779  apop_data * apop_model_hessian(apop_data * data, apop_model *model, double delta) ;
780 #else
781  apop_data * apop_model_hessian_base(apop_data * data, apop_model *model, double delta) ;
782  apop_varad_declare(apop_data *, apop_model_hessian, apop_data * data; apop_model *model; double delta);
783 #define apop_model_hessian(...) apop_varad_link(apop_model_hessian, __VA_ARGS__)
784 #endif
785 
786 #ifdef APOP_NO_VARIADIC
787  apop_data * apop_model_numerical_covariance(apop_data * data, apop_model *model, double delta) ;
788 #else
789  apop_data * apop_model_numerical_covariance_base(apop_data * data, apop_model *model, double delta) ;
790  apop_varad_declare(apop_data *, apop_model_numerical_covariance, apop_data * data; apop_model *model; double delta);
791 #define apop_model_numerical_covariance(...) apop_varad_link(apop_model_numerical_covariance, __VA_ARGS__)
792 #endif
793 
794 
795 void apop_maximum_likelihood(apop_data * data, apop_model *dist);
796 
797 #ifdef APOP_NO_VARIADIC
798  apop_model * apop_estimate_restart (apop_model *e, apop_model *copy, char * starting_pt, double boundary) ;
799 #else
800  apop_model * apop_estimate_restart_base(apop_model *e, apop_model *copy, char * starting_pt, double boundary) ;
801  apop_varad_declare(apop_model *, apop_estimate_restart, apop_model *e; apop_model *copy; char * starting_pt; double boundary);
802 #define apop_estimate_restart(...) apop_varad_link(apop_estimate_restart, __VA_ARGS__)
803 #endif
804 
805 
806 //in apop_linear_constraint.c
807 #ifdef APOP_NO_VARIADIC
808  long double apop_linear_constraint(gsl_vector *beta, apop_data * constraint, double margin) ;
809 #else
810  long double apop_linear_constraint_base(gsl_vector *beta, apop_data * constraint, double margin) ;
811  apop_varad_declare(long double, apop_linear_constraint, gsl_vector *beta; apop_data * constraint; double margin);
812 #define apop_linear_constraint(...) apop_varad_link(apop_linear_constraint, __VA_ARGS__)
813 #endif
814 
815 
816 //in apop_model_fix_params.c
819 
820 
821 
823 
825 /* This declares the vtable macros for each procedure that uses the mechanism.
826 
827 --We want to have type-checking on the functions put into the vtables. Type checking
828 happens only with functions, not macros, so we need a type_check function for every
829 vtable.
830 
831 --Only once in your codebase, you'll need to #define Declare_type_checking_fns to
832 actually define the type checking function. Everywhere else, the function is merely
833 declared.
834 
835 --All other uses point to having a macro, such as using __VA_ARGS__ to allow any sort
836 of inputs to the hash.
837 
838 --We want to have such a macro for every vtable. That means that we need a macro
839 to write macros. We can't do that with C macros, so this file uses m4 macros to
840 generate C macros.
841 
842 --After the m4 definition of make_vtab_fns, each new vtable requires a typedef, a hash
843 definition, and a call to make_vtab_fns to do the rest.
844 */
845 
846 
847 int apop_vtable_add(char const *tabname, void *fn_in, unsigned long hash);
848 void *apop_vtable_get(char const *tabname, unsigned long hash);
849 int apop_vtable_drop(char const *tabname, unsigned long hash);
850 
851 typedef apop_model *(*apop_update_type)(apop_data *, apop_model* , apop_model*);
852 #define apop_update_hash(m1, m2) ( \
853  ((m1)->log_likelihood ? (size_t)(m1)->log_likelihood : \
854  (m1)->p ? (size_t)(m1)->p*33 : \
855  (m1)->draw ? (size_t)(m1)->draw*33*27 \
856  : 33*27*19) \
857  +((m2)->log_likelihood ? (size_t)(m2)->log_likelihood : \
858  (m2)->p ? (size_t)(m2)->p*33 : \
859  (m2)->draw ? (size_t)(m2)->draw*33*27 \
860  : 33*27*19 \
861  ) * 37)
862 #ifdef Declare_type_checking_fns
863 void apop_update_type_check(apop_update_type in){ };
864 #else
865 void apop_update_type_check(apop_update_type in);
866 #endif
867 #define apop_update_vtable_add(fn, ...) apop_update_type_check(fn), apop_vtable_add("apop_update", fn, apop_update_hash(__VA_ARGS__))
868 #define apop_update_vtable_get(...) apop_vtable_get("apop_update", apop_update_hash(__VA_ARGS__))
869 #define apop_update_vtable_drop(...) apop_vtable_drop("apop_update", apop_update_hash(__VA_ARGS__))
870 
871 typedef long double (*apop_entropy_type)(apop_model *model);
872 #define apop_entropy_hash(m1) ((size_t)(m1)->log_likelihood + 33 * (size_t)((m1)->p) + 27*(size_t)((m1)->draw))
873 #ifdef Declare_type_checking_fns
874 void apop_entropy_type_check(apop_entropy_type in){ };
875 #else
876 void apop_entropy_type_check(apop_entropy_type in);
877 #endif
878 #define apop_entropy_vtable_add(fn, ...) apop_entropy_type_check(fn), apop_vtable_add("apop_entropy", fn, apop_entropy_hash(__VA_ARGS__))
879 #define apop_entropy_vtable_get(...) apop_vtable_get("apop_entropy", apop_entropy_hash(__VA_ARGS__))
880 #define apop_entropy_vtable_drop(...) apop_vtable_drop("apop_entropy", apop_entropy_hash(__VA_ARGS__))
881 
882 typedef void (*apop_score_type)(apop_data *d, gsl_vector *gradient, apop_model *params);
883 #define apop_score_hash(m1) ((size_t)((m1)->log_likelihood ? (m1)->log_likelihood : (m1)->p))
884 #ifdef Declare_type_checking_fns
885 void apop_score_type_check(apop_score_type in){ };
886 #else
887 void apop_score_type_check(apop_score_type in);
888 #endif
889 #define apop_score_vtable_add(fn, ...) apop_score_type_check(fn), apop_vtable_add("apop_score", fn, apop_score_hash(__VA_ARGS__))
890 #define apop_score_vtable_get(...) apop_vtable_get("apop_score", apop_score_hash(__VA_ARGS__))
891 #define apop_score_vtable_drop(...) apop_vtable_drop("apop_score", apop_score_hash(__VA_ARGS__))
892 
893 typedef apop_model* (*apop_parameter_model_type)(apop_data *, apop_model *);
894 #define apop_parameter_model_hash(m1) ((size_t)((m1)->log_likelihood ? (m1)->log_likelihood : (m1)->p)*33 + (m1)->estimate ? (size_t)(m1)->estimate: 27)
895 #ifdef Declare_type_checking_fns
896 void apop_parameter_model_type_check(apop_parameter_model_type in){ };
897 #else
898 void apop_parameter_model_type_check(apop_parameter_model_type in);
899 #endif
900 #define apop_parameter_model_vtable_add(fn, ...) apop_parameter_model_type_check(fn), apop_vtable_add("apop_parameter_model", fn, apop_parameter_model_hash(__VA_ARGS__))
901 #define apop_parameter_model_vtable_get(...) apop_vtable_get("apop_parameter_model", apop_parameter_model_hash(__VA_ARGS__))
902 #define apop_parameter_model_vtable_drop(...) apop_vtable_drop("apop_parameter_model", apop_parameter_model_hash(__VA_ARGS__))
903 
904 typedef apop_data * (*apop_predict_type)(apop_data *d, apop_model *params);
905 #define apop_predict_hash(m1) ((size_t)((m1)->log_likelihood ? (m1)->log_likelihood : (m1)->p)*33 + (m1)->estimate ? (size_t)(m1)->estimate: 27)
906 #ifdef Declare_type_checking_fns
907 void apop_predict_type_check(apop_predict_type in){ };
908 #else
909 void apop_predict_type_check(apop_predict_type in);
910 #endif
911 #define apop_predict_vtable_add(fn, ...) apop_predict_type_check(fn), apop_vtable_add("apop_predict", fn, apop_predict_hash(__VA_ARGS__))
912 #define apop_predict_vtable_get(...) apop_vtable_get("apop_predict", apop_predict_hash(__VA_ARGS__))
913 #define apop_predict_vtable_drop(...) apop_vtable_drop("apop_predict", apop_predict_hash(__VA_ARGS__))
914 
915 typedef void (*apop_model_print_type)(apop_model *params, FILE *out);
916 #define apop_model_print_hash(m1) ((m1)->log_likelihood ? (size_t)(m1)->log_likelihood : \
917  (m1)->p ? (size_t)(m1)->p*33 : \
918  (m1)->estimate ? (size_t)(m1)->estimate*33*33 : \
919  (m1)->draw ? (size_t)(m1)->draw*33*27 : \
920  (m1)->cdf ? (size_t)(m1)->cdf*27*27 \
921  : 27)
922 #ifdef Declare_type_checking_fns
923 void apop_model_print_type_check(apop_model_print_type in){ };
924 #else
925 void apop_model_print_type_check(apop_model_print_type in);
926 #endif
927 #define apop_model_print_vtable_add(fn, ...) apop_model_print_type_check(fn), apop_vtable_add("apop_model_print", fn, apop_model_print_hash(__VA_ARGS__))
928 #define apop_model_print_vtable_get(...) apop_vtable_get("apop_model_print", apop_model_print_hash(__VA_ARGS__))
929 #define apop_model_print_vtable_drop(...) apop_vtable_drop("apop_model_print", apop_model_print_hash(__VA_ARGS__))
930  //End of Doxygen ignore.
932 
933 
935 
936 long double apop_generalized_harmonic(int N, double s) __attribute__ ((__pure__));
937 
938 apop_data * apop_test_anova_independence(apop_data *d);
939 #define apop_test_ANOVA_independence(d) apop_test_anova_independence(d)
940 
941 #ifdef APOP_NO_VARIADIC
942  int apop_regex(const char *string, const char* regex, apop_data **substrings, const char use_case) ;
943 #else
944  int apop_regex_base(const char *string, const char* regex, apop_data **substrings, const char use_case) ;
945  apop_varad_declare(int, apop_regex, const char *string; const char* regex; apop_data **substrings; const char use_case);
946 #define apop_regex(...) apop_varad_link(apop_regex, __VA_ARGS__)
947 #endif
948 
949 
950 int apop_system(const char *fmt, ...) __attribute__ ((format (printf,1,2)));
951 
952 //Histograms and PMFs
953 gsl_vector * apop_vector_moving_average(gsl_vector *, size_t);
954 apop_data * apop_histograms_test_goodness_of_fit(apop_model *h0, apop_model *h1);
955 apop_data * apop_test_kolmogorov(apop_model *m1, apop_model *m2);
956 apop_data *apop_data_pmf_compress(apop_data *in);
957 #ifdef APOP_NO_VARIADIC
958  apop_data * apop_data_to_bins(apop_data const *indata, apop_data const *binspec, int bin_count, char close_top_bin) ;
959 #else
960  apop_data * apop_data_to_bins_base(apop_data const *indata, apop_data const *binspec, int bin_count, char close_top_bin) ;
961  apop_varad_declare(apop_data *, apop_data_to_bins, apop_data const *indata; apop_data const *binspec; int bin_count; char close_top_bin);
962 #define apop_data_to_bins(...) apop_varad_link(apop_data_to_bins, __VA_ARGS__)
963 #endif
964 
965 #ifdef APOP_NO_VARIADIC
966  apop_model * apop_model_to_pmf(apop_model *model, apop_data *binspec, long int draws, int bin_count) ;
967 #else
968  apop_model * apop_model_to_pmf_base(apop_model *model, apop_data *binspec, long int draws, int bin_count) ;
969  apop_varad_declare(apop_model *, apop_model_to_pmf, apop_model *model; apop_data *binspec; long int draws; int bin_count);
970 #define apop_model_to_pmf(...) apop_varad_link(apop_model_to_pmf, __VA_ARGS__)
971 #endif
972 
973 
974 //text conveniences
975 #ifdef APOP_NO_VARIADIC
976  char* apop_text_paste(apop_data const*strings, char *between, char *before, char *after, char *between_cols, int (*prune)(apop_data* , int , int , void*), void* prune_parameter) ;
977 #else
978  char* apop_text_paste_base(apop_data const*strings, char *between, char *before, char *after, char *between_cols, int (*prune)(apop_data* , int , int , void*), void* prune_parameter) ;
979  apop_varad_declare(char*, apop_text_paste, apop_data const*strings; char *between; char *before; char *after; char *between_cols; int (*prune)(apop_data* , int , int , void*); void* prune_parameter);
980 #define apop_text_paste(...) apop_varad_link(apop_text_paste, __VA_ARGS__)
981 #endif
982 
993 #define Apop_notify(verbosity, ...) {\
994  if (apop_opts.verbose != -1 && apop_opts.verbose >= verbosity) { \
995  if (!apop_opts.log_file) apop_opts.log_file = stderr; \
996  fprintf(apop_opts.log_file, "%s: ", __func__); fprintf(apop_opts.log_file, __VA_ARGS__); fprintf(apop_opts.log_file, "\n"); \
997  fflush(apop_opts.log_file); \
998 } }
999 
1001 #define Apop_maybe_abort(level) \
1002  {if ((apop_opts.verbose >= level && apop_opts.stop_on_warning == 'v') \
1003  || (apop_opts.stop_on_warning=='w') ) \
1004  raise(SIGTRAP);}
1005 
1030 #define Apop_stopif(test, onfail, level, ...) do {\
1031  if (test) { \
1032  Apop_notify(level, __VA_ARGS__); \
1033  Apop_maybe_abort(level) \
1034  onfail; \
1035  } } while(0)
1036 
1037 #define apop_errorlevel -5
1038 
1040 //For use in stopif, to return a blank apop_data set with an error attached.
1041 #define apop_return_data_error(E) {apop_data *out=apop_data_alloc(); out->error='E'; return out;}
1042 
1043 /* The Apop_stopif macro is currently favored, but there's a long history of prior
1044  error-handling setups. Consider all of the Assert... macros below to be deprecated.
1045 */
1046 #define Apop_assert_c(test, returnval, level, ...) \
1047  Apop_stopif(!(test), return returnval, level, __VA_ARGS__)
1048 
1049 #define Apop_assert(test, ...) Apop_assert_c((test), 0, apop_errorlevel, __VA_ARGS__)
1050 
1051 //For things that return void. Transitional and deprecated at birth.
1052 #define Apop_assert_n(test, ...) Apop_assert_c((test), , apop_errorlevel, __VA_ARGS__)
1053 #define Apop_assert_negone(test, ...) Apop_assert_c((test), -1, apop_errorlevel, __VA_ARGS__)
1054  //End of Doxygen ignore.
1055 
1056 //Missing data
1057 #ifdef APOP_NO_VARIADIC
1058  apop_data * apop_data_listwise_delete(apop_data *d, char inplace) ;
1059 #else
1060  apop_data * apop_data_listwise_delete_base(apop_data *d, char inplace) ;
1061  apop_varad_declare(apop_data *, apop_data_listwise_delete, apop_data *d; char inplace);
1062 #define apop_data_listwise_delete(...) apop_varad_link(apop_data_listwise_delete, __VA_ARGS__)
1063 #endif
1064 
1065 apop_model * apop_ml_impute(apop_data *d, apop_model* meanvar);
1066 
1067 #ifdef APOP_NO_VARIADIC
1068  apop_model *apop_model_metropolis(apop_data *d, gsl_rng* rng, apop_model *m);
1069 #else
1070  apop_model * apop_model_metropolis_base(apop_data *d, gsl_rng* rng, apop_model *m);
1071  apop_varad_declare(apop_model *, apop_model_metropolis, apop_data *d; gsl_rng* rng; apop_model *m);
1072 #define apop_model_metropolis(...) apop_varad_link(apop_model_metropolis, __VA_ARGS__)
1073 #endif
1074 
1075 #ifdef APOP_NO_VARIADIC
1076  apop_model * apop_update(apop_data *data, apop_model *prior, apop_model *likelihood, gsl_rng *rng) ;
1077 #else
1078  apop_model * apop_update_base(apop_data *data, apop_model *prior, apop_model *likelihood, gsl_rng *rng) ;
1079  apop_varad_declare(apop_model *, apop_update, apop_data *data; apop_model *prior; apop_model *likelihood; gsl_rng *rng);
1080 #define apop_update(...) apop_varad_link(apop_update, __VA_ARGS__)
1081 #endif
1082 
1083 
1084 #ifdef APOP_NO_VARIADIC
1085  double apop_test(double statistic, char *distribution, double p1, double p2, char tail) ;
1086 #else
1087  double apop_test_base(double statistic, char *distribution, double p1, double p2, char tail) ;
1088  apop_varad_declare(double, apop_test, double statistic; char *distribution; double p1; double p2; char tail);
1089 #define apop_test(...) apop_varad_link(apop_test, __VA_ARGS__)
1090 #endif
1091 
1092 
1093 //apop_sort.c
1094 #ifdef APOP_NO_VARIADIC
1095  apop_data *apop_data_sort(apop_data *data, apop_data *sort_order, char asc, char inplace, double *col_order);
1096 #else
1097  apop_data * apop_data_sort_base(apop_data *data, apop_data *sort_order, char asc, char inplace, double *col_order);
1098  apop_varad_declare(apop_data *, apop_data_sort, apop_data *data; apop_data *sort_order; char asc; char inplace; double *col_order);
1099 #define apop_data_sort(...) apop_varad_link(apop_data_sort, __VA_ARGS__)
1100 #endif
1101 
1102 
1103 //raking
1104 #ifdef APOP_NO_VARIADIC
1105  apop_data * apop_rake(char const *margin_table, char * const*var_list,
1106  int var_ct, char * const *contrasts, int contrast_ct,
1107  char const *structural_zeros, int max_iterations, double tolerance,
1108  char const *count_col, char const *init_table,
1109  char const *init_count_col, double nudge) ;
1110 #else
1111  apop_data * apop_rake_base(char const *margin_table, char * const*var_list,
1112  int var_ct, char * const *contrasts, int contrast_ct,
1113  char const *structural_zeros, int max_iterations, double tolerance,
1114  char const *count_col, char const *init_table,
1115  char const *init_count_col, double nudge) ;
1116  apop_varad_declare(apop_data *, apop_rake, char const *margin_table; char * const*var_list;
1117  int var_ct; char * const *contrasts; int contrast_ct;
1118  char const *structural_zeros; int max_iterations; double tolerance;
1119  char const *count_col; char const *init_table;
1120  char const *init_count_col; double nudge);
1121 #define apop_rake(...) apop_varad_link(apop_rake, __VA_ARGS__)
1122 #endif
1123 
1124 
1125 
1126 #include <gsl/gsl_cdf.h>
1127 #include <gsl/gsl_blas.h>
1128 #include <gsl/gsl_sf_log.h>
1129 #include <gsl/gsl_sf_exp.h>
1130 #include <gsl/gsl_linalg.h>
1131 #include <gsl/gsl_sf_gamma.h>
1132 #include <gsl/gsl_sf_psi.h>
1133 #include <gsl/gsl_randist.h>
1134 #include <gsl/gsl_histogram.h>
1135 #include <gsl/gsl_statistics_double.h>
1136 
1137 
1138  //Some linear algebra utilities
1139 
1140 double apop_det_and_inv(const gsl_matrix *in, gsl_matrix **out, int calc_det, int calc_inv);
1141 #ifdef APOP_NO_VARIADIC
1142  apop_data * apop_dot(const apop_data *d1, const apop_data *d2, char form1, char form2) ;
1143 #else
1144  apop_data * apop_dot_base(const apop_data *d1, const apop_data *d2, char form1, char form2) ;
1145  apop_varad_declare(apop_data *, apop_dot, const apop_data *d1; const apop_data *d2; char form1; char form2);
1146 #define apop_dot(...) apop_varad_link(apop_dot, __VA_ARGS__)
1147 #endif
1148 
1149 #ifdef APOP_NO_VARIADIC
1150  int apop_vector_bounded(const gsl_vector *in, long double max) ;
1151 #else
1152  int apop_vector_bounded_base(const gsl_vector *in, long double max) ;
1153  apop_varad_declare(int, apop_vector_bounded, const gsl_vector *in; long double max);
1154 #define apop_vector_bounded(...) apop_varad_link(apop_vector_bounded, __VA_ARGS__)
1155 #endif
1156 
1157 gsl_matrix * apop_matrix_inverse(const gsl_matrix *in) ;
1158 double apop_matrix_determinant(const gsl_matrix *in) ;
1159 //apop_data* apop_sv_decomposition(gsl_matrix *data, int dimensions_we_want);
1160 #ifdef APOP_NO_VARIADIC
1161  apop_data * apop_matrix_pca(gsl_matrix *data, int const dimensions_we_want) ;
1162 #else
1163  apop_data * apop_matrix_pca_base(gsl_matrix *data, int const dimensions_we_want) ;
1164  apop_varad_declare(apop_data *, apop_matrix_pca, gsl_matrix *data; int const dimensions_we_want);
1165 #define apop_matrix_pca(...) apop_varad_link(apop_matrix_pca, __VA_ARGS__)
1166 #endif
1167 
1168 #ifdef APOP_NO_VARIADIC
1169  gsl_vector * apop_vector_stack(gsl_vector *v1, gsl_vector const * v2, char inplace) ;
1170 #else
1171  gsl_vector * apop_vector_stack_base(gsl_vector *v1, gsl_vector const * v2, char inplace) ;
1172  apop_varad_declare(gsl_vector *, apop_vector_stack, gsl_vector *v1; gsl_vector const * v2; char inplace);
1173 #define apop_vector_stack(...) apop_varad_link(apop_vector_stack, __VA_ARGS__)
1174 #endif
1175 
1176 #ifdef APOP_NO_VARIADIC
1177  gsl_matrix * apop_matrix_stack(gsl_matrix *m1, gsl_matrix const * m2, char posn, char inplace) ;
1178 #else
1179  gsl_matrix * apop_matrix_stack_base(gsl_matrix *m1, gsl_matrix const * m2, char posn, char inplace) ;
1180  apop_varad_declare(gsl_matrix *, apop_matrix_stack, gsl_matrix *m1; gsl_matrix const * m2; char posn; char inplace);
1181 #define apop_matrix_stack(...) apop_varad_link(apop_matrix_stack, __VA_ARGS__)
1182 #endif
1183 
1184 
1185 void apop_vector_log(gsl_vector *v);
1186 void apop_vector_log10(gsl_vector *v);
1187 void apop_vector_exp(gsl_vector *v);
1188 
1190 
1193 #define APOP_SUBMATRIX(m, srow, scol, nrows, ncols, o) gsl_matrix apop_mm_##o = gsl_matrix_submatrix((m), (srow), (scol), (nrows),(ncols)).matrix;\
1194 gsl_matrix * o = &( apop_mm_##o ); // Use \ref Apop_subm.
1195 #define Apop_submatrix APOP_SUBMATRIX
1196 
1197 #define Apop_col_v(m, col, v) gsl_vector apop_vv_##v = ((col) == -1) ? (gsl_vector){} : gsl_matrix_column((m)->matrix, (col)).vector;\
1198 gsl_vector * v = ((col)==-1) ? (m)->vector : &( apop_vv_##v ); // Use \ref Apop_cv.
1199 
1200 #define Apop_row_v(m, row, v) Apop_matrix_row((m)->matrix, row, v) // Use \ref Apop_rv.
1201 #define Apop_rows(d, rownum, len, outd) apop_data *outd = Apop_rs(d, rownum, len) // Use \ref Apop_rs.
1202 #define Apop_row(d, row, outd) Apop_rows(d, row, 1, outd) // Use \ref Apop_r.
1203 #define Apop_cols(d, colnum, len, outd) apop_data *outd = Apop_cs(d, colnum, len); // Use \ref Apop_cs.
1204  //End of Doxygen ignore.
1205 
1211 #define Apop_row_tv(m, row, v) gsl_vector apop_vv_##v = gsl_matrix_row((m)->matrix, apop_name_find((m)->names, row, 'r')).vector;\
1212 gsl_vector * v = &( apop_vv_##v );
1213 
1220 #define Apop_col_tv(m, col, v) gsl_vector apop_vv_##v = gsl_matrix_column((m)->matrix, apop_name_find((m)->names, col, 'c')).vector;\
1221 gsl_vector * v = &( apop_vv_##v );
1222 
1228 #define Apop_row_t(d, rowname, outd) int apop_row_##outd = apop_name_find((d)->names, rowname, 'r'); Apop_rows(d, apop_row_##outd, 1, outd)
1229 
1235 #define Apop_col_t(d, colname, outd) int apop_col_##outd = apop_name_find((d)->names, colname, 'c'); Apop_cols(d, apop_col_##outd, 1, outd)
1236 
1237 // The above versions relied on gsl_views, which stick to C as of 1989 CE.
1238 // Better to just create the views via designated initializers.
1239 
1240 
1251 #define Apop_subm(matrix_to_view, srow, scol, nrows, ncols)( \
1252  (!(matrix_to_view) \
1253  || (matrix_to_view)->size1 < (srow)+(nrows) || (srow) < 0 \
1254  || (matrix_to_view)->size2 < (scol)+(ncols) || (scol) < 0) ? NULL \
1255  : &(gsl_matrix){.size1=(nrows), .size2=(ncols), \
1256  .tda=(matrix_to_view)->tda, \
1257  .data=gsl_matrix_ptr((matrix_to_view), (srow), (scol))} \
1258  )
1259 
1270 #define Apop_mrv(matrix_to_view, row) Apop_rv(&(apop_data){.matrix=matrix_to_view}, row)
1271 
1287 #define Apop_mcv(matrix_to_view, col) Apop_cv(&(apop_data){.matrix=matrix_to_view}, col)
1288 
1303 #define Apop_rv(data_to_view, row) ( \
1304  ((data_to_view) == NULL || (data_to_view)->matrix == NULL \
1305  || (data_to_view)->matrix->size1 <= (row) || (row) < 0) ? NULL \
1306  : &(gsl_vector){.size=(data_to_view)->matrix->size2, \
1307  .stride=1, .data=gsl_matrix_ptr((data_to_view)->matrix, (row), 0)} \
1308  )
1309 
1328 #define Apop_cv(data_to_view, col) ( \
1329  !(data_to_view) ? NULL \
1330  : (col)==-1 ? (data_to_view)->vector \
1331  : (!(data_to_view)->matrix \
1332  || (data_to_view)->matrix->size2 <= (col) || ((int)(col)) < -1) ? NULL \
1333  : &(gsl_vector){.size=(data_to_view)->matrix->size1, \
1334  .stride=(data_to_view)->matrix->tda, .data=gsl_matrix_ptr((data_to_view)->matrix, 0, (col))} \
1335  )
1336 
1338 /* Not (yet) for public use. */
1339 #define Apop_subvector(v, start, len) ( \
1340  ((v) == NULL || (v)->size < ((start)+(len)) || (start) < 0) ? NULL \
1341  : &(gsl_vector){.size=(len), .stride=(v)->stride, .data=(v)->data+(start*(v)->stride)})
1342 
1351 #define Apop_rs(d, rownum, len)( \
1352  (!(d) || (rownum) < 0) ? NULL \
1353  : &(apop_data){ \
1354  .names= ( !((d)->names) ? NULL : \
1355  &(apop_name){ \
1356  .title = (d)->names->title, \
1357  .vector = (d)->names->vector, \
1358  .col = (d)->names->col, \
1359  .row = ((d)->names->row && (d)->names->rowct > (rownum)) ? &((d)->names->row[rownum]) : NULL, \
1360  .text = (d)->names->text, \
1361  .colct = (d)->names->colct, \
1362  .rowct = (d)->names->row ? (GSL_MIN(1, GSL_MAX((d)->names->rowct - (int)(rownum), 0))) \
1363  : 0, \
1364  .textct = (d)->names->textct }), \
1365  .vector= Apop_subvector((d->vector), (rownum), (len)), \
1366  .matrix = Apop_subm(((d)->matrix), (rownum), 0, (len), (d)->matrix?(d)->matrix->size2:0), \
1367  .weights = Apop_subvector(((d)->weights), (rownum), (len)), \
1368  .textsize[0]=(d)->textsize[0]> (rownum)+(len)-1 ? (len) : 0, \
1369  .textsize[1]=(d)->textsize[1], \
1370  .text = (d)->text ? &((d)->text[rownum]) : NULL, \
1371  })
1372 
1373 
1380 #define Apop_cs(d, colnum, len) ( \
1381  (!(d)||!(d)->matrix || (d)->matrix->size2 <= (colnum)+(len)-1) \
1382  ? NULL \
1383  : &(apop_data){ \
1384  .vector= NULL, \
1385  .weights= (d)->weights, \
1386  .matrix = Apop_subm((d)->matrix, 0, colnum, (d)->matrix->size1, (len)),\
1387  .textsize[0] = 0, \
1388  .textsize[1] = 0, \
1389  .text = NULL, \
1390  .names= (d)->names ? &(apop_name){ \
1391  .title = (d)->names->title, \
1392  .vector = NULL, \
1393  .row = (d)->names->row, \
1394  .col = ((d)->names->col && (d)->names->colct > colnum) ? &((d)->names->col[colnum]) : NULL, \
1395  .text = NULL, \
1396  .rowct = (d)->names->rowct, \
1397  .colct = (d)->names->col ? (GSL_MIN(len, GSL_MAX((d)->names->colct - colnum, 0))) \
1398  : 0, \
1399  .textct = (d)->names->textct } : NULL \
1400  })
1401 
1420 #define Apop_r(d, rownum) Apop_rs(d, rownum, 1)
1421 
1429 #define Apop_c(d, col) Apop_cs(d, col, 1)
1430 
1432 #define APOP_COL Apop_col
1433 #define apop_col Apop_col
1434 #define APOP_COL_T Apop_col_t
1435 #define apop_col_t Apop_col_t
1436 #define APOP_COL_TV Apop_col_tv
1437 #define apop_col_tv Apop_col_tv
1438 
1439 #define APOP_ROW Apop_row
1440 #define apop_row Apop_row
1441 #define APOP_COLS Apop_cols
1442 #define apop_cols Apop_cols
1443 #define APOP_COL_V Apop_col_v
1444 #define apop_col_v Apop_col_v
1445 #define APOP_ROW_V Apop_row_v
1446 #define apop_row_v Apop_row_v
1447 #define APOP_ROWS Apop_rows
1448 #define apop_rows Apop_rows
1449 #define Apop_data_row Apop_row #deprecated
1450 #define APOP_ROW_T Apop_row_t
1451 #define apop_row_t Apop_row_t
1452 #define APOP_ROW_TV Apop_row_tv
1453 #define apop_row_tv Apop_row_tv
1454 
1456 #define Apop_matrix_row(m, row, v) gsl_vector apop_vv_##v = gsl_matrix_row((m), (row)).vector;\
1457 gsl_vector * v = &( apop_vv_##v );
1458 
1459 /* Deprecated. Use Apop_mcv */
1460 #define Apop_matrix_col(m, col, v) gsl_vector apop_vv_##v = gsl_matrix_column((m), (col)).vector;\
1461 gsl_vector * v = &( apop_vv_##v );
1462 
1463 #define APOP_MATRIX_ROW Apop_matrix_row
1464 #define apop_matrix_row Apop_matrix_row
1465 #define APOP_MATRIX_COL Apop_matrix_col
1466 #define apop_matrix_col Apop_matrix_col
1467 
1470 long double apop_vector_sum(const gsl_vector *in);
1471 double apop_vector_var_m(const gsl_vector *in, const double mean);
1472 #ifdef APOP_NO_VARIADIC
1473  double apop_vector_correlation(const gsl_vector *ina, const gsl_vector *inb, const gsl_vector *weights) ;
1474 #else
1475  double apop_vector_correlation_base(const gsl_vector *ina, const gsl_vector *inb, const gsl_vector *weights) ;
1476  apop_varad_declare(double, apop_vector_correlation, const gsl_vector *ina; const gsl_vector *inb; const gsl_vector *weights);
1477 #define apop_vector_correlation(...) apop_varad_link(apop_vector_correlation, __VA_ARGS__)
1478 #endif
1479 
1480 double apop_vector_kurtosis(const gsl_vector *in);
1481 double apop_vector_skew(const gsl_vector *in);
1482 
1483 #define apop_sum apop_vector_sum
1484 #define apop_var apop_vector_var
1485 #define apop_mean apop_vector_mean
1486 
1488 
1489 #ifdef APOP_NO_VARIADIC
1490  int apop_table_exists(char const *name, char remove) ;
1491 #else
1492  int apop_table_exists_base(char const *name, char remove) ;
1493  apop_varad_declare(int, apop_table_exists, char const *name; char remove);
1494 #define apop_table_exists(...) apop_varad_link(apop_table_exists, __VA_ARGS__)
1495 #endif
1496 
1497 
1498 int apop_db_open(char const *filename);
1499 #ifdef APOP_NO_VARIADIC
1500  int apop_db_close(char vacuum) ;
1501 #else
1502  int apop_db_close_base(char vacuum) ;
1503  apop_varad_declare(int, apop_db_close, char vacuum);
1504 #define apop_db_close(...) apop_varad_link(apop_db_close, __VA_ARGS__)
1505 #endif
1506 
1507 
1508 int apop_query(const char *q, ...) __attribute__ ((format (printf,1,2)));
1509 apop_data * apop_query_to_text(const char * fmt, ...) __attribute__ ((format (printf,1,2)));
1510 apop_data * apop_query_to_data(const char * fmt, ...) __attribute__ ((format (printf,1,2)));
1511 apop_data * apop_query_to_mixed_data(const char *typelist, const char * fmt, ...) __attribute__ ((format (printf,2,3)));
1512 gsl_vector * apop_query_to_vector(const char * fmt, ...) __attribute__ ((format (printf,1,2)));
1513 double apop_query_to_float(const char * fmt, ...) __attribute__ ((format (printf,1,2)));
1514 
1515 int apop_data_to_db(const apop_data *set, const char *tabname, char);
1516 
1517 
1519 
1520  //Part I: macros and fns for getting/setting settings groups and elements
1521 
1523 void * apop_settings_get_grp(apop_model *m, char *type, char fail);
1524 void apop_settings_remove_group(apop_model *m, char *delme);
1525 void apop_settings_copy_group(apop_model *outm, apop_model *inm, char *copyme);
1526 void *apop_settings_group_alloc(apop_model *model, char *type, void *free_fn, void *copy_fn, void *the_group);
1527 apop_model *apop_settings_group_alloc_wm(apop_model *model, char *type, void *free_fn, void *copy_fn, void *the_group); //End of Doxygen ignore.
1529 
1543 #define Apop_settings_get_group(m, type) apop_settings_get_grp(m, #type, 'c')
1544 
1549 #define Apop_settings_rm_group(m, type) apop_settings_remove_group(m, #type)
1550 
1563 #define Apop_settings_add_group(model, type, ...) \
1564  apop_settings_group_alloc(model, #type, type ## _settings_free, type ## _settings_copy, type ##_settings_init ((type ## _settings) {__VA_ARGS__}))
1565 
1570 #define apop_model_copy_set(model, type, ...) \
1571  apop_settings_group_alloc_wm(apop_model_copy(model), #type, type ## _settings_free, type ## _settings_copy, type ##_settings_init ((type ## _settings) {__VA_ARGS__}))
1572 
1573 
1589 #define Apop_model_set_settings(model, ...) \
1590  apop_settings_group_alloc_wm(apop_model_copy(model), #model, model ## _settings_free, model ## _settings_copy, model ##_settings_init ((model ## _settings) {__VA_ARGS__}))
1591 
1592 #define apop_model_set_settings Apop_model_set_settings
1593 
1600 #define Apop_settings_get(model, type, setting) \
1601  (((type ## _settings *) apop_settings_get_grp(model, #type, 'f'))->setting)
1602 
1608 #define Apop_settings_set(model, type, setting, data) \
1609  do { \
1610  if (!(model)) continue; /* silent fail. */ \
1611  type ## _settings *apop_tmp_settings = apop_settings_get_grp(model, #type, 'c'); \
1612  Apop_stopif(!apop_tmp_settings, (model)->error='s', 0, "You're trying to modify a setting in " \
1613  #model "'s setting group of type " #type " but that model doesn't have such a group."); \
1614  apop_tmp_settings->setting = (data); \
1615  } while (0);
1616 
1618 #define Apop_settings_add Apop_settings_set
1619 #define APOP_SETTINGS_ADD Apop_settings_set
1620 #define apop_settings_set Apop_settings_set
1621 #define APOP_SETTINGS_GET Apop_settings_get
1622 #define apop_settings_get Apop_settings_get
1623 #define APOP_SETTINGS_ADD_GROUP Apop_settings_add_group
1624 #define apop_settings_add_group Apop_settings_add_group
1625 #define APOP_SETTINGS_GET_GROUP Apop_settings_get_group
1626 #define apop_settings_get_group Apop_settings_get_group
1627 #define APOP_SETTINGS_RM_GROUP Apop_settings_rm_group
1628 #define apop_settings_rm_group Apop_settings_rm_group
1629 #define Apop_model_copy_set apop_model_copy_set
1630 
1631 //deprecated:
1632 #define Apop_model_add_group Apop_settings_add_group
1633  //End of Doxygen ignore.
1635 
1639 #define Apop_settings_declarations(ysg) \
1640  ysg##_settings * ysg##_settings_init(ysg##_settings); \
1641  void * ysg##_settings_copy(ysg##_settings *); \
1642  void ysg##_settings_free(ysg##_settings *);
1643 
1647 #define Apop_settings_init(name, ...) \
1648  name##_settings *name##_settings_init(name##_settings in) { \
1649  name##_settings *out = malloc(sizeof(name##_settings)); \
1650  *out = in; \
1651  __VA_ARGS__; \
1652  return out; \
1653  }
1654 
1656 #define Apop_varad_set(var, value) (out)->var = (in).var ? (in).var : (value);
1657 
1662 #define Apop_settings_copy(name, ...) \
1663  void * name##_settings_copy(name##_settings *in) {\
1664  name##_settings *out = malloc(sizeof(name##_settings)); \
1665  *out = *in; \
1666  __VA_ARGS__; \
1667  return out; \
1668  }
1669 
1673 #define Apop_settings_free(name, ...) \
1674  void name##_settings_free(name##_settings *in) {\
1675  __VA_ARGS__; \
1676  free(in); \
1677  }
1678 
1679  //Part II: the details of extant settings groups.
1680 
1681 
1683 typedef struct{
1684  double *starting_pt;
1687  char *method;
1710  double step_size,
1711  tolerance,
1712 delta;
1713  int max_iterations;
1715  int verbose;
1717  double dim_cycle_tolerance;
1724 //simulated annealing (also uses step_size);
1725  int n_tries, iters_fixed_T;
1726  double k, t_initial, mu_t, t_min ;
1727  gsl_rng *rng;
1728  apop_data **path;
1733 
1735 typedef struct {
1736  int destroy_data;
1737  apop_data *instruments;
1738  char want_cov;
1739  char want_expected_value;
1740  apop_model *input_distribution;
1742 
1766 typedef struct {
1767  //init/copy/free are in apop_mle.c
1768  char covariance; /*< If 'y', calculate the covariance matrix. Default 'n'. */
1769  char predicted;/*< If 'y', calculate the predicted values. This is typically as many
1770  items as rows in your data set. Default 'n'. */
1771  char tests;/*< If 'y', run any hypothesis tests offered by the model's estimation routine. Default 'n'. */
1772  char info;/*< If 'y', add an info table with elements such as log likelihood or AIC. Default 'n'. */
1774 
1776 typedef struct {
1777  int draws;
1778  gsl_rng *rng;
1779  gsl_matrix *draws_made;
1780  int *draws_refcount;
1782 
1783 
1785 typedef struct {
1786  apop_model *base;
1787  int index;
1788  gsl_rng *rng;
1789  int draws;
1791 
1792 
1794 typedef struct {
1795  gsl_vector *cmf;
1796  char draw_index;
1798  long double total_weight;
1799  int *cmf_refct;
1801 
1802 
1804 typedef struct{
1805  apop_data *base_data;
1806  apop_model *base_pmf;
1808  apop_model *kernel;
1810  void (*set_fn)(apop_data*, apop_model*);
1815  int own_pmf, own_kernel;
1817 
1818 struct apop_mcmc_settings;
1819 
1833 typedef struct apop_mcmc_proposal_s {
1834  apop_model *proposal;
1841  void (*step_fn)(double const *, struct apop_mcmc_proposal_s*, struct apop_mcmc_settings *);
1849  int (*adapt_fn)(struct apop_mcmc_proposal_s *ps, struct apop_mcmc_settings *ms);
1853  int accept_count, reject_count;
1858 
1860 typedef struct apop_mcmc_settings {
1861  apop_data *data;
1862  long int periods;
1863  double burnin;
1865  int histosegments;
1866  double last_ll;
1867  apop_model *pmf;
1869  apop_model *base_model;
1872  apop_mcmc_proposal_s *proposals;
1875  int proposal_count;
1876  double target_accept_rate;
1877  int accept_count;
1878  int reject_count;
1879  char gibbs_chunks;
1892  size_t *block_starts;
1893  int block_count, proposal_is_cp;
1895  char start_at;
1900  void (*base_step_fn)(double const *, struct apop_mcmc_proposal_s*, struct apop_mcmc_settings *);
1901  int (*base_adapt_fn)(struct apop_mcmc_proposal_s *ps, struct apop_mcmc_settings *ms);
1904 
1906 //Loess, including the old FORTRAN-to-C.
1907 struct loess_struct {
1908  struct {
1909  long n, p;
1910  double *y, *x;
1911  double *weights;
1912  } in;
1913  struct {
1914  double span;
1915  long degree;
1916  long normalize;
1917  long parametric[8];
1918  long drop_square[8];
1919  char *family;
1920  } model;
1921  struct {
1922  char *surface;
1923  char *statistics;
1924  double cell;
1925  char *trace_hat;
1926  long iterations;
1927  } control;
1928  struct {
1929  long *parameter, *a;
1930  double *xi, *vert, *vval;
1931  } kd_tree;
1932  struct {
1933  double *fitted_values;
1934  double *fitted_residuals;
1935  double enp, s;
1936  double one_delta, two_delta;
1937  double *pseudovalues;
1938  double trace_hat;
1939  double *diagonal;
1940  double *robust;
1941  double *divisor;
1942  } out;
1943 }; //End of Doxygen ignore.
1945 
1958 typedef struct {
1959  apop_data *data;
1960  struct loess_struct lo_s;
2049  int want_predict_ci;
2051  double ci_level;
2054 
2055 
2057 typedef struct point { /* a point in the x,y plane */
2058  double x,y; /* x and y coordinates */
2059  double ey; /* exp(y-ymax+YCEIL) */
2060  double cum; /* integral up to x of rejection envelope */
2061  int f; /* is y an evaluated point of log-density */
2062  struct point *pl,*pr; /* envelope points to left and right of x */
2063 } POINT;
2064 
2065 /* This includes the envelope info and the metropolis steps. */
2066 typedef struct { /* attributes of the entire rejection envelope */
2067  int cpoint; /* number of POINTs in current envelope */
2068  int npoint; /* max number of POINTs allowed in envelope */
2069  double ymax; /* the maximum y-value in the current envelope */
2070  POINT *p; /* start of storage of envelope POINTs */
2071  double *convex; /* adjustment for convexity */
2072  double metro_xprev; /* previous Markov chain iterate */
2073  double metro_yprev; /* current log density at xprev */
2074 } arms_state;
2085 typedef struct {
2086  double *xinit;
2090  double xl;
2091  double xr;
2092  double convex;
2093  int ninit;
2094  int npoint;
2095  char do_metro;
2097  double xprev;
2098  int neval;
2099  arms_state *state;
2100  apop_model *model;
2102 
2103 
2105 typedef struct {
2106  char *splitpage;
2107  apop_model *model1;
2108  apop_model *model2;
2110 
2111 typedef struct {
2112  apop_data *(*base_to_transformed)(apop_data*);
2113  apop_data *(*transformed_to_base)(apop_data*);
2114  double (*jacobian_to_base)(apop_data*);
2115  apop_model *base_model;
2121 typedef struct {
2122  apop_model *base_model;
2123  double (*constraint)(apop_data *, apop_model *);
2124  double (*scaling)(apop_model *);
2125  gsl_rng *rng;
2127  int draw_ct;
2129 
2130 typedef struct {
2131  apop_model *generator_m;
2132  apop_model *ll_m;
2133  int draw_ct;
2140 typedef struct {
2141  gsl_vector *weights;
2142  apop_model **model_list;
2143  int model_count;
2144  int *param_sizes;
2145  apop_model *cmf;
2146  int *cmf_refct;
2148 
2149  //Models built via call to apop_model_copy_set.
2150 
2151 #define apop_model_dcompose(...) Apop_model_set_settings(apop_composition, __VA_ARGS__)
2152 #define apop_model_dconstrain(...) Apop_model_set_settings(apop_dconstrain, __VA_ARGS__)
2153 #define apop_model_coordinate_transform(...) Apop_model_set_settings(apop_coordinate_transform, __VA_ARGS__)
2154 
2155 //Doxygen drops whatever is after these declarations, so I put them last.
2161 Apop_settings_declarations(apop_arms)
2162 Apop_settings_declarations(apop_mcmc)
2163 Apop_settings_declarations(apop_loess)
2164 Apop_settings_declarations(apop_cross)
2165 Apop_settings_declarations(apop_mixture)
2166 Apop_settings_declarations(apop_dconstrain)
2167 Apop_settings_declarations(apop_composition)
2168 Apop_settings_declarations(apop_parts_wanted)
2169 Apop_settings_declarations(apop_kernel_density)
2170 Apop_settings_declarations(apop_coordinate_transform)
2171 
2172 #ifdef __cplusplus
2173 }
2174 #endif
2175  //End doxygen's all_public grouping
2177 
2178 //Part of the intent of a convenience header like this is that you
2179 //don't have to remember what else you're including. So here are
2180 //some other common GSL headers:
2181 #include <math.h>
2182 #include <gsl/gsl_sort.h>
2183 #include <gsl/gsl_eigen.h>
2184 #include <gsl/gsl_sort_vector.h>
2185 #include <gsl/gsl_permutation.h>
2186 #include <gsl/gsl_integration.h>
apop_data * apop_data_transpose(apop_data *in, char transpose_text, char inplace)
Definition: apop_data.c:1152
void apop_text_free(char ***freeme, int rows, int cols)
Definition: apop_data.c:142
apop_data * apop_estimate_coefficient_of_determination(apop_model *)
Definition: apop_regression.c:544
apop_model * apop_estimate_restart(apop_model *e, apop_model *copy, char *starting_pt, double boundary)
Definition: apop_mle.c:731
Definition: docs/include/apop.h:72
apop_opts_type apop_opts
Definition: docs/include/apop.h:142
apop_data * apop_data_rm_page(apop_data *data, const char *title, const char free_p)
Definition: apop_data.c:1412
apop_data * apop_data_add_page(apop_data *dataset, apop_data *newpage, const char *title)
Definition: apop_data.c:1371
apop_data * apop_data_alloc(const size_t size1, const size_t size2, const int size3)
Definition: apop_data.c:34
apop_model * apop_model_clear(apop_data *data, apop_model *model)
Definition: apop_model.c:21
apop_data * apop_data_get_page(const apop_data *data, const char *title, const char match)
Definition: apop_data.c:1338
apop_model * apop_ml_impute(apop_data *d, apop_model *meanvar)
Definition: apop_missing_data.c:140
apop_data * apop_rake(char const *margin_table, char *const *var_list, int var_ct, char *const *contrasts, int contrast_ct, char const *structural_zeros, int max_iterations, double tolerance, char const *count_col, char const *init_table, char const *init_count_col, double nudge)
Definition: apop_rake.c:527
long double apop_generalized_harmonic(int N, double s)
Definition: apop_asst.c:146
int apop_regex(const char *string, const char *regex, apop_data **substrings, const char use_case)
Definition: apop_asst.c:256
int apop_name_find(const apop_name *n, const char *findme, const char type)
Definition: apop_name.c:203
gsl_vector * apop_matrix_map(const gsl_matrix *m, double(*fn)(gsl_vector *))
Definition: apop_mapply.c:337
apop_model * apop_beta_from_mean_var(double m, double v)
Definition: apop_asst.c:338
void apop_data_add_named_elmt(apop_data *d, char *name, double val)
Definition: apop_data.c:984
apop_dconstrain
Definition: model_doc.h:396
apop_data * apop_histograms_test_goodness_of_fit(apop_model *h0, apop_model *h1)
Definition: apop_hist.c:70
int apop_text_to_db(char const *text_file, char *tabname, int has_row_names, int has_col_names, char **field_names, int const *field_ends, apop_data *field_params, char *table_params, char const *delimiters, char if_table_exists)
Definition: apop_conversions.c:1168
gsl_vector * apop_numerical_gradient(apop_data *data, apop_model *model, double delta)
Definition: apop_mle.c:128
void apop_vector_log10(gsl_vector *v)
Definition: apop_linear_algebra.c:162
apop_data * apop_data_rank_compress(apop_data *in, int min_bins)
Definition: apop_conversions.c:283
apop_binomial
Definition: model_doc.h:865
gsl_vector * apop_vector_map(const gsl_vector *v, double(*fn)(double))
Definition: apop_mapply.c:373
double apop_vector_map_sum(const gsl_vector *in, double(*fn)(double))
Definition: apop_mapply.c:443
apop_data * apop_model_hessian(apop_data *data, apop_model *model, double delta)
Definition: apop_mle.c:196
void apop_score(apop_data *d, gsl_vector *out, apop_model *m)
Definition: apop_model.c:298
int apop_arms_draw(double *out, gsl_rng *r, apop_model *m)
Definition: apop_arms.c:89
apop_cross
Definition: model_doc.h:732
double * apop_vector_percentiles(gsl_vector *data, char rounding)
Definition: apop_stats.c:486
apop_beta
Definition: model_doc.h:62
double apop_matrix_to_positive_semidefinite(gsl_matrix *m)
Definition: apop_stats.c:966
apop_data * apop_test_fisher_exact(apop_data *intab)
Definition: apop_fexact.c:1894
apop_mixture
Definition: model_doc.h:470
apop_data * apop_text_alloc(apop_data *in, const size_t row, const size_t col)
Definition: apop_data.c:1063
apop_loess
Definition: model_doc.h:119
apop_data * apop_model_draws(apop_model *model, int count, apop_data *draws)
Definition: apop_asst.c:421
Definition: docs/include/apop.h:1794
void apop_vector_print(gsl_vector *data, Output_declares)
Definition: apop_output.c:207
apop_data * apop_data_calloc(const size_t size1, const size_t size2, const int size3)
Definition: apop_data.c:86
apop_data * apop_data_correlation(const apop_data *in)
Definition: apop_stats.c:666
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)
Definition: apop_mapply.c:494
apop_data * apop_data_to_factors(apop_data *data, char intype, int incol, int outcol)
Definition: apop_regression.c:439
double apop_vector_skew_pop(gsl_vector const *v, gsl_vector const *weights)
Definition: apop_stats.c:90
apop_name * apop_name_copy(apop_name *in)
Definition: apop_name.c:184
double apop_log_likelihood(apop_data *d, apop_model *m)
Definition: apop_model.c:265
apop_name * apop_name_alloc(void)
Definition: apop_name.c:24
long double apop_multivariate_lngamma(double a, int p)
Definition: apop_stats.c:887
double * apop_data_ptr(apop_data *data, int row, int col, const char *rowname, const char *colname, const char *page)
Definition: apop_data.c:786
int apop_table_exists(char const *name, char remove)
Definition: apop_db.c:132
gsl_matrix * apop_matrix_stack(gsl_matrix *m1, gsl_matrix const *m2, char posn, char inplace)
Definition: apop_linear_algebra.c:258
gsl_matrix * apop_vector_to_matrix(const gsl_vector *in, char row_col)
Definition: apop_conversions.c:70
apop_logit
Definition: model_doc.h:975
apop_data * apop_data_stack(apop_data *m1, apop_data *m2, char posn, char inplace)
Definition: apop_data.c:365
Definition: docs/include/apop.h:118
gsl_vector * apop_vector_stack(gsl_vector *v1, gsl_vector const *v2, char inplace)
Definition: apop_linear_algebra.c:195
Definition: docs/include/apop.h:1735
apop_data * apop_data_copy(const apop_data *in)
Definition: apop_data.c:295
void apop_crosstab_to_db(apop_data *in, char *tabname, char *row_col_name, char *col_col_name, char *data_col_name)
Definition: apop_conversions.c:213
double apop_test(double statistic, char *distribution, double p1, double p2, char tail)
Definition: apop_tests.c:448
void apop_vector_log(gsl_vector *v)
Definition: apop_linear_algebra.c:170
double apop_cdf(apop_data *d, apop_model *m)
Definition: apop_model.c:535
Definition: docs/include/apop.h:1860
apop_data * apop_dot(const apop_data *d1, const apop_data *d2, char form1, char form2)
Definition: apop_linear_algebra.c:410
void apop_prep(apop_data *d, apop_model *m)
Definition: apop_model.c:455
apop_data * apop_data_rm_rows(apop_data *in, int *drop, int(*do_drop)(apop_data *, void *), void *drop_parameter)
Definition: apop_data.c:1475
long double apop_linear_constraint(gsl_vector *beta, apop_data *constraint, double margin)
Definition: apop_linear_constraint.c:126
apop_coordinate_transform
Definition: model_doc.h:695
apop_data * apop_matrix_pca(gsl_matrix *data, int const dimensions_we_want)
Definition: apop_linear_algebra.c:112
double apop_vector_cov(gsl_vector const *v1, gsl_vector const *v2, gsl_vector const *weights)
Definition: apop_stats.c:594
apop_bernoulli
Definition: model_doc.h:7
double apop_matrix_mean(const gsl_matrix *data)
Definition: apop_stats.c:360
apop_t_distribution
Definition: model_doc.h:214
void apop_name_free(apop_name *free_me)
Definition: apop_name.c:115
apop_data * apop_db_to_crosstab(char const *tabname, char const *row, char const *col, char const *data, char is_aggregate)
Definition: apop_conversions.c:126
Definition: docs/include/apop.h:85
apop_data * apop_bootstrap_cov(apop_data *data, apop_model *model, gsl_rng *rng, int iterations, char keep_boots, char ignore_nans, apop_data **boot_store)
Definition: apop_bootstrap.c:128
apop_probit
Definition: model_doc.h:634
apop_data * apop_data_to_bins(apop_data const *indata, apop_data const *binspec, int bin_count, char close_top_bin)
Definition: apop_hist.c:333
apop_model * apop_model_metropolis(apop_data *d, gsl_rng *rng, apop_model *m)
Definition: apop_mcmc.c:283
int(* base_adapt_fn)(struct apop_mcmc_proposal_s *ps, struct apop_mcmc_settings *ms)
Definition: docs/include/apop.h:1901
apop_model * apop_parameter_model(apop_data *d, apop_model *m)
Definition: apop_model.c:364
int apop_db_close(char vacuum)
Definition: apop_db.c:181
double apop_vector_kurtosis(const gsl_vector *in)
Definition: apop_stats.c:53
Definition: docs/include/apop.h:2130
Definition: docs/include/apop.h:1958
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)
Definition: apop_mapply.c:121
int apop_name_add(apop_name *n, char const *add_me, char type)
Definition: apop_name.c:42
char apop_data_free_base(apop_data *freeme)
Definition: apop_data.c:166
Definition: docs/include/apop.h:1766
void apop_vector_normalize(gsl_vector *in, gsl_vector **out, const char normalization_type)
Definition: apop_stats.c:298
void apop_vector_exp(gsl_vector *v)
Definition: apop_linear_algebra.c:178
double apop_vector_correlation(const gsl_vector *ina, const gsl_vector *inb, const gsl_vector *weights)
Definition: apop_stats.c:178
Definition: docs/include/apop.h:98
Definition: docs/include/apop.h:2140
apop_data * apop_text_unique_elements(const apop_data *d, size_t col)
Definition: apop_regression.c:90
apop_data * apop_test_kolmogorov(apop_model *m1, apop_model *m2)
Definition: apop_hist.c:249
double apop_vector_kurtosis_pop(gsl_vector const *v, gsl_vector const *weights)
Definition: apop_stats.c:129
void apop_model_free(apop_model *free_me)
Definition: apop_model.c:56
void apop_matrix_mean_and_var(const gsl_matrix *data, double *mean, double *var)
Definition: apop_stats.c:391
gsl_vector * apop_vector_realloc(gsl_vector *v, size_t newheight)
Definition: apop_data.c:1309
double apop_matrix_map_sum(const gsl_matrix *in, double(*fn)(gsl_vector *))
Definition: apop_mapply.c:471
apop_data * apop_data_pmf_compress(apop_data *in)
Definition: apop_pmf.c:340
long double apop_multivariate_gamma(double a, int p)
Definition: apop_stats.c:875
apop_improper_uniform
Definition: model_doc.h:610
void apop_model_print(apop_model *model, FILE *output_pipe)
Definition: apop_model.c:121
apop_data * apop_text_to_data(char const *text_file, int has_row_names, int has_col_names, int const *field_ends, char const *delimiters)
Definition: apop_conversions.c:629
double apop_vector_distance(const gsl_vector *ina, const gsl_vector *inb, const char metric, const double norm)
Definition: apop_stats.c:221
apop_model * apop_update(apop_data *data, apop_model *prior, apop_model *likelihood, gsl_rng *rng)
Definition: apop_update.c:185
apop_kernel_density
Definition: model_doc.h:902
Definition: docs/include/apop.h:2121
apop_data * apop_anova(char *table, char *data, char *grouping1, char *grouping2)
Definition: apop_tests.c:317
int apop_text_set(apop_data *in, const size_t row, const size_t col, const char *fmt,...)
Definition: apop_data.c:1035
apop_multivariate_normal
Definition: model_doc.h:80
Definition: docs/include/apop.h:62
apop_model * apop_model_copy(apop_model *in)
Definition: apop_model.c:163
Definition: docs/include/apop.h:1833
gsl_matrix * apop_matrix_inverse(const gsl_matrix *in)
Definition: apop_linear_algebra.c:71
void apop_data_rm_columns(apop_data *d, int *drop)
Definition: apop_data.c:696
apop_model * apop_model_to_pmf(apop_model *model, apop_data *binspec, long int draws, int bin_count)
Definition: apop_hist.c:31
apop_data * apop_data_prune_columns_base(apop_data *d, char **colnames)
Definition: apop_data.c:738
apop_exponential
Definition: model_doc.h:671
Definition: docs/include/apop.h:2111
Definition: docs/include/apop.h:2105
apop_normal
Definition: model_doc.h:242
double apop_vector_mean(gsl_vector const *v, gsl_vector const *weights)
Definition: apop_stats.c:521
apop_data * apop_data_covariance(const apop_data *in)
Definition: apop_stats.c:642
apop_dirichlet
Definition: model_doc.h:266
int apop_matrix_is_positive_semidefinite(gsl_matrix *m, char semi)
Definition: apop_stats.c:934
apop_data * apop_jackknife_cov(apop_data *data, apop_model *model)
Definition: apop_bootstrap.c:51
int apop_vector_bounded(const gsl_vector *in, long double max)
Definition: apop_linear_algebra.c:331
apop_data ** apop_data_split(apop_data *in, int splitpoint, char r_or_c)
Definition: apop_data.c:478
int apop_draw(double *out, gsl_rng *r, apop_model *m)
Definition: apop_model.c:421
void apop_vector_apply(gsl_vector *v, void(*fn)(double *))
Definition: apop_mapply.c:390
#define Apop_settings_declarations(ysg)
Definition: apop.h:1639
Definition: docs/include/apop.h:1785
int accept_count
Definition: docs/include/apop.h:1877
apop_zipf
Definition: model_doc.h:37
long double apop_vector_entropy(gsl_vector *in)
Definition: apop_stats.c:703
apop_uniform
Definition: model_doc.h:713
apop_lognormal
Definition: model_doc.h:162
apop_data * apop_data_sort(apop_data *data, apop_data *sort_order, char asc, char inplace, double *col_order)
Definition: apop_sort.c:144
apop_model * apop_model_fix_params_get_base(apop_model *model_in)
Definition: apop_fix_params.c:275
Definition: docs/include/apop.h:2085
gsl_matrix * apop_matrix_map_all(const gsl_matrix *in, double(*fn)(double))
Definition: apop_mapply.c:408
gsl_vector * apop_array_to_vector(double *in, int size)
Definition: apop_conversions.c:39
void apop_name_print(apop_name *n)
Definition: apop_name.c:84
void apop_matrix_apply(gsl_matrix *m, void(*fn)(gsl_vector *))
Definition: apop_mapply.c:355
apop_data * apop_paired_t_test(gsl_vector *a, gsl_vector *b)
Definition: apop_tests.c:92
apop_model * apop_model_fix_params(apop_model *model_in)
Definition: apop_fix_params.c:232
gsl_matrix * apop_matrix_copy(const gsl_matrix *in)
Definition: apop_conversions.c:356
double apop_data_get(const apop_data *data, size_t row, int col, const char *rowname, const char *colname, const char *page)
Definition: apop_data.c:841
apop_model * apop_estimate(apop_data *d, apop_model *m)
Definition: apop_model.c:237
apop_pmf
Definition: model_doc.h:324
double apop_vector_skew(const gsl_vector *in)
Definition: apop_stats.c:41
apop_yule
Definition: model_doc.h:949
gsl_rng * apop_rng_alloc(int seed)
Definition: apop_bootstrap.c:19
apop_data * apop_model_numerical_covariance(apop_data *data, apop_model *model, double delta)
Definition: apop_mle.c:255
apop_iv
Definition: model_doc.h:844
int reject_count
Definition: docs/include/apop.h:1878
gsl_matrix * apop_matrix_realloc(gsl_matrix *m, size_t newheight, size_t newwidth)
Definition: apop_data.c:1259
void apop_data_memcpy(apop_data *out, const apop_data *in)
Definition: apop_data.c:216
int apop_db_open(char const *filename)
Definition: apop_db.c:86
apop_data * apop_predict(apop_data *d, apop_model *m)
Definition: apop_model.c:490
long double apop_model_entropy(apop_model *in, int draws)
Definition: apop_stats.c:748
double apop_matrix_map_all_sum(const gsl_matrix *in, double(*fn)(double))
Definition: apop_mapply.c:457
apop_gamma
Definition: model_doc.h:807
apop_poisson
Definition: model_doc.h:1028
apop_multinomial
Definition: model_doc.h:753
long double apop_matrix_sum(const gsl_matrix *m)
Definition: apop_stats.c:346
apop_data * apop_test_anova_independence(apop_data *d)
Definition: apop_tests.c:243
double apop_p(apop_data *d, apop_model *m)
Definition: apop_model.c:250
long double apop_kl_divergence(apop_model *from, apop_model *to, int draw_ct, gsl_rng *rng)
Definition: apop_stats.c:810
void apop_maximum_likelihood(apop_data *data, apop_model *dist)
Definition: apop_mle.c:657
gsl_vector * apop_vector_copy(const gsl_vector *in)
Definition: apop_conversions.c:338
double apop_rng_GHgB3(gsl_rng *r, double *a)
Definition: apop_asst.c:313
Definition: docs/include/apop.h:1804
double apop_matrix_determinant(const gsl_matrix *in)
Definition: apop_linear_algebra.c:85
double apop_det_and_inv(const gsl_matrix *in, gsl_matrix **out, int calc_det, int calc_inv)
Definition: apop_linear_algebra.c:41
Definition: docs/include/apop.h:1683
gsl_vector * apop_data_pack(const apop_data *in, gsl_vector *out, char more_pages, char use_info_pages)
Definition: apop_conversions.c:805
Definition: docs/include/apop.h:1776
apop_data * apop_data_to_dummies(apop_data *d, int col, char type, int keep_first, char append, char remove)
Definition: apop_regression.c:335
apop_data * apop_data_rank_expand(apop_data *in)
Definition: apop_conversions.c:314
apop_ols
Definition: model_doc.h:510
apop_data * apop_data_get_factor_names(apop_data *data, int col, char type)
Definition: apop_regression.c:155
char * apop_text_paste(apop_data const *strings, char *between, char *before, char *after, char *between_cols, int(*prune)(apop_data *, int, int, void *), void *prune_parameter)
Definition: apop_asst.c:100
void apop_data_print(const apop_data *data, Output_declares)
Definition: apop_output.c:333
double apop_vector_var_m(const gsl_vector *in, const double mean)
Definition: apop_stats.c:159
void apop_matrix_apply_all(gsl_matrix *in, void(*fn)(double *))
Definition: apop_mapply.c:428
void apop_matrix_print(const gsl_matrix *data, Output_declares)
Definition: apop_output.c:365
long double apop_vector_sum(const gsl_vector *in)
Definition: apop_stats.c:16
int apop_data_set(apop_data *data, size_t row, int col, const double val, const char *rowname, const char *colname, const char *page)
Definition: apop_data.c:915
void apop_name_stack(apop_name *n1, apop_name *nadd, char type1, char typeadd)
Definition: apop_name.c:139
void apop_data_unpack(const gsl_vector *in, apop_data *d, char use_info_pages)
Definition: apop_conversions.c:727
apop_data * apop_data_summarize(apop_data *data)
Definition: apop_stats.c:422
apop_data * apop_t_test(gsl_vector *a, gsl_vector *b)
Definition: apop_tests.c:55
double apop_vector_var(gsl_vector const *v, gsl_vector const *weights)
Definition: apop_stats.c:559
gsl_vector * apop_vector_unique_elements(const gsl_vector *v)
Definition: apop_regression.c:63
apop_data * apop_f_test(apop_model *est, apop_data *contrast)
Definition: apop_tests.c:140
apop_data * apop_data_listwise_delete(apop_data *d, char inplace)
Definition: apop_missing_data.c:35
gsl_vector * apop_vector_moving_average(gsl_vector *, size_t)
Definition: apop_hist.c:382
Autogenerated by doxygen (Debian ).