Patterns in static

Apophenia

Some examples

Here are a few pieces of sample code, gathered from elsewhere in the documentation, for testing your installation or to give you a sense of what code with Apophenia's tools looks like. If you'd like more context or explanation, please click through to the page from which the example was taken.

In the documentation for the apop_ols model, a program to read in data and run a regression. You'll need to go to that page for the sample data and further discussion.

From the Setting up page, an example of gathering data from two processes, saving the input to a database, then doing a later analysis:

#include <apop.h>
//Your processes are probably a bit more complex.
double process_one(gsl_rng *r){
return gsl_rng_uniform(r) * gsl_rng_uniform(r) ;
}
double process_two(gsl_rng *r){
return gsl_rng_uniform(r);
}
int main(){
double p1, p2;
int i;
gsl_rng *r = apop_rng_alloc(123);
//create the database and the data table.
apop_db_open("runs.db");
apop_table_exists("samples", 'd'); //If the table already exists, delete it.
apop_query("create table samples(iteration, process, value); begin;");
//populate the data table with runs.
for (i=0; i<1000; i++){
p1 = process_one(r);
p2 = process_two(r);
apop_query("insert into samples values(%i, %i, %g);", i, 1, p1);
apop_query("insert into samples values(%i, %i, %g);", i, 2, p2);
}
apop_query("commit;"); //the begin-commit wrapper saves writes to the drive.
//pull the data from the database, converting it into a table along the way.
apop_data *m = apop_db_to_crosstab("samples", "iteration","process", "value");
gsl_vector *v1 = Apop_cv(m, 0); //get vector views of the two table columns.
gsl_vector *v2 = Apop_cv(m, 1);
//Output a table of means/variances, and t-test results.
printf("\t mean\t\t var\n");
printf("process 1: %f\t%f\n", apop_mean(v1), apop_var(v1));
printf("process 2: %f\t%f\n\n", apop_mean(v2), apop_var(v2));
printf("t test\n");
apop_data_print(m, "the_data.txt");
}

In the An outline of the library section on map/apply, a new $t$-test on every row, with all operations acting on entire rows rather than individual data points:

#include <apop.h>
double row_offset;
gsl_rng *r;
void offset_rng(double *v){*v = gsl_rng_uniform(r) + row_offset;}
double find_tstat(gsl_vector *in){ return apop_mean(in)/sqrt(apop_var(in));}
double conf(double in, void *df){ return gsl_cdf_tdist_P(in, *(int *)df);}
//apop_vector_mean is a macro, so we can't point a pointer to it.
double mu(gsl_vector *in){ return apop_vector_mean(in);}
int main(){
apop_data *d = apop_data_alloc(10, 100);
r = apop_rng_alloc(3242);
for (int i=0; i< 10; i++){
row_offset = gsl_rng_uniform(r)*2 -1;
Apop_row_v(d, i, onerow);
apop_vector_apply(onerow, offset_rng);
}
size_t df = d->matrix->size2-1;
apop_data *means = apop_map(d, .fn_v = mu, .part ='r');
apop_data *tstats = apop_map(d, .fn_v = find_tstat, .part ='r');
apop_data *confidences = apop_map(tstats, .fn_dp = conf, .param = &df);
#ifndef Testing
printf("means:\n"); apop_data_show(means);
printf("\nt stats:\n"); apop_data_show(tstats);
printf("\nconfidences:\n"); apop_data_show(confidences);
#endif
//Some sanity checks, for Apophenia's test suite.
for (int i=0; i< 10; i++){
//sign of mean == sign of t stat.
assert(apop_data_get(means, i, -1) * apop_data_get(tstats, i, -1) >=0);
//inverse of P-value should be the t statistic.
assert(fabs(gsl_cdf_tdist_Pinv(apop_data_get(confidences, i, -1),100)
- apop_data_get(tstats, i, -1) < 1e-3));
}
}

In the documentation for apop_query_to_text, a program to list all the tables in an SQLite database.

#include <apop.h>
void print_table_list(char *db_file){
apop_db_open(db_file);
apop_data *tab_list= apop_query_to_text("select name "
"from sqlite_master where type=='table'");
for(int i=0; i< tab_list->textsize[0]; i++)
printf("%s\n", tab_list->text[i][0]);
}
int main(int argc, char **argv){
if (argc == 1){
printf("Give me a database name, and I will print out "
"the list of tables contained therein.\n");
return 0;
}
print_table_list(argv[1]);
}

A demonstration of fixing parameters to create a marginal distribution, via apop_model_fix_params

#include <apop.h>
int main(){
size_t ct = 5e4;
//set up the model & params
apop_data *params = apop_data_alloc(2,2,2);
apop_data_fill(params, 8, 1, 0.5,
2, 0.5, 1);
pvm->parameters = apop_data_copy(params);
pvm->dsize = 2;
apop_data *d = apop_model_draws(pvm, ct);
//set up and estimate a model with fixed covariance matrix but free means
gsl_vector_set_all(pvm->parameters->vector, GSL_NAN);
apop_model *e1 = apop_estimate(d, mep1);
//compare results, via assert for the test suite, or on-screen for human use.
#ifdef Testing
assert(apop_vector_distance(params->vector, e1->parameters->vector)<1e-2);
#else
printf("original params: ");
apop_vector_show(params->vector);
printf("estimated params: ");
#endif
}

Several uses of the apop_dot function

/* A demonstration of dot products and various useful
transformations among types. */
#include <apop.h>
double eps=1e-3;//slow to converge series-->large tolerance.
#define Diff(L, R) Apop_assert(fabs((L)-(R)<(eps)), "%g is too different from %g (abitrary limit=%g).", (double)(L), (double)(R), eps);
int main(){
int len = 3000;
gsl_vector *v = gsl_vector_alloc(len);
for (double i=0; i< len; i++) gsl_vector_set(v, i, 1./(i+1));
double square;
gsl_blas_ddot(v, v, &square);
#ifndef Testing
printf("1 + (1/2)^2 + (1/3)^2 + ...= %g\n", square);
#endif
double pi_over_six = gsl_pow_2(M_PI)/6.;
Diff(square, pi_over_six);
/* Now using apop_dot, in a few forms.
First, vector-as-data dot itself.
If one of the inputs is a vector,
apop_dot puts the output in a vector-as-data:*/
apop_data *v_as_data = apop_vector_to_data(v);
apop_data *vdotv = apop_dot(v_as_data, v_as_data);
Diff(gsl_vector_get(vdotv->vector, 0), pi_over_six);
/* As a matrix-as-data, via a compound literal.
Default is a len X 1 matrix. */
gsl_matrix *v_as_matrix = apop_vector_to_matrix(v);
apop_data dm = (apop_data){.matrix=v_as_matrix};
// (1 X len) dot (len X 1) --- OK, produce a scalar.
apop_data *mdotv = apop_dot(v_as_data, &dm);
double scalarval = gsl_vector_get(mdotv->vector, 0);
Diff(scalarval, pi_over_six);
//(len X 1) dot (len X 1) --- bad dimensions.
apop_opts.verbose=-1; //don't print an error.
apop_data *mdotv2 = apop_dot(&dm, v_as_data);
apop_opts.verbose=0; //back to sanity.
assert(mdotv2->error);
// If we want (len X 1) dot (1 X len) --> (len X len),
// use apop_vector_to_matrix.
apop_data dmr = (apop_data){.matrix=apop_vector_to_matrix(v, .row_col='r')};
apop_data *product_matrix = apop_dot(&dm, &dmr);
//The trace is the sum of squares:
gsl_vector_view trace = gsl_matrix_diagonal(product_matrix->matrix);
double tracesum = apop_sum(&trace.vector);
Diff(tracesum, pi_over_six);
apop_data_free(product_matrix);
gsl_matrix_free(dmr.matrix);
}

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