![]() |
|
Macros | |
#define | omp_threadnum 0 |
Typedefs | |
typedef int(* | apop_fn_riip )(apop_data *, int, int, void *) |
Functions | |
char * | apop_text_paste (apop_data const *strings, char *between, char *before, char *after, char *between_cols, apop_fn_riip prune, void *prune_parameter) |
double | apop_generalized_harmonic (int N, double s) |
int | apop_system (const char *fmt,...) |
double * | apop_vector_percentiles (gsl_vector *data, char rounding) |
int | apop_regex (const char *string, const char *regex, apop_data **substrings, const char use_case) |
double | apop_rng_GHgB3 (gsl_rng *r, double *a) |
apop_model * | apop_beta_from_mean_var (double m, double v) |
gsl_rng * | apop_rng_get_thread_base (int thread) |
apop_data * | apop_model_draws (apop_model *model, int count, apop_data *draws) |
Variables | |
char * | apop_nul_string |
The odds and ends bin. Copyright (c) 2005–2007, 2010 by Ben Klemens. Licensed under the modified GNU GPL v2; see COPYING and COPYING2.
apop_model* apop_beta_from_mean_var | ( | double | m, |
double | v | ||
) |
The Beta distribution is useful for modeling because it is bounded between zero and one, and can be either unimodal (if the variance is low) or bimodal (if the variance is high), and can have either a slant toward the bottom or top of the range (depending on the mean).
The distribution has two parameters, typically named and
, which can be difficult to interpret. However, there is a one-to-one mapping between (alpha, beta) pairs and (mean, variance) pairs. Since we have good intuition about the meaning of means and variances, this function takes in a mean and variance, calculates alpha and beta behind the scenes, and returns a random draw from the appropriate Beta distribution.
m | The mean the Beta distribution should have. Notice that m is in [0,1]. |
v | The variance which the Beta distribution should have. It is in (0, 1/12), where (1/12) is the variance of a Uniform(0,1) distribution. Funny things happen with variance near 1/12 and mean far from 1/2. |
apop_beta
model with its parameters appropriately set. out->error=='r' | Range error: mean is not within [0, 1]. |
double apop_generalized_harmonic | ( | int | N, |
double | s | ||
) |
Calculate
N
is zero or negative, return NaN. Notify the user if apop_opts.verbosity >=1
For example:
apop_data* apop_model_draws | ( | apop_model * | model, |
int | count, | ||
apop_data * | draws | ||
) |
Make a set of random draws from a model and write them to an apop_data set.
model | The model from which draws will be made. Must already be prepared and/or estimated. |
count | The number of draws to make. If draw_matrix is not NULL , then this is ignored and count=draw_matrix->matrix->size1 . default=1000. |
draws | If not NULL , a pre-allocated data set whose matrix element will be filled with draws. |
size
draws. If draw_matrix!=NULL
, then return a pointer to it.out->error=='m' | Input model isn't good for making draws: it is NULL , or m->dsize=0 . |
out->error=='s' | You gave me a draws matrix, but its size is less than the size of a single draw from the data, model->dsize . |
out->error=='d' | Trouble drawing from the distribution for at least one row. That row is set to all NAN . |
NULL apop_data
set, but its matrix
element is NULL
, when apop_opts.verbose>=1
.Here is a two-line program to draw a different set of ten Standard Normals on every run (provided runs are more than a second apart):
int apop_regex | ( | const char * | string, |
const char * | regex, | ||
apop_data ** | substrings, | ||
const char | use_case | ||
) |
A convenience function for regular expression searching
For example, "p.val" will match "P value", "p.value", "p values" (and even "tempeval", so be careful).
If you give a non-NULL
address in which to place a table of paren-delimited substrings, I'll return them as a row in the text element of the returned apop_data set. I'll return all the matches, filling the first row with substrings from the first application of your regex, then filling the next row with another set of matches (if any), and so on to the end of the string. Useful when parsing a list of items, for example.
string | The string to search (no default) |
regex | The regular expression (no default) |
substrings | Parens in the regex indicate that I should return matching substrings. Give me the address of an apop_data* set, and I will allocate and fill the text portion with matches. Default= NULL , meaning do not return substrings (even if parens exist in the regex). If no match, return an empty apop_data set, so output->textsize[0]==0 . |
use_case | Should I be case sensitive, 'y' or 'n' ? (default = 'n' , which is not the POSIX default.) |
substrings
may be allocated and filled if needed.apop_opts.stop_on_warning='n'
returns -1 on error (e.g., regex NULL
or didn't compile). strings==NULL
, I return 0—no match—and if substrings
is provided, set it to NULL
.&subs
, not plain subs
. Also, the non-match has a zero-length blank in subs->text[0][1]
.([A-Za-z])([0-9])
, the column zero of outdata
will hold letters, and column one will hold numbers. Use apop_data_transpose to reverse this so that the letters are in outdata->text[0]
and numbers in outdata->text[1]
. double apop_rng_GHgB3 | ( | gsl_rng * | r, |
double * | a | ||
) |
RNG from a Generalized Hypergeometric type B3.
Devroye uses this as the base for many of his distribution-generators, including the Waring.
GSL_NAN
if the function doesn't stop. int apop_system | ( | const char * | fmt, |
... | |||
) |
char* apop_text_paste | ( | apop_data const * | strings, |
char * | between, | ||
char * | before, | ||
char * | after, | ||
char * | between_cols, | ||
apop_fn_riip | prune, | ||
void * | prune_parameter | ||
) |
Join together a list or array of strings, with optional separators between the strings.
strings | An apop_data set with a grid of text to be combined into a single string |
between | The text to put in between the rows of the table, such as ", ". (Default is a single space: " ") |
before | The text to put at the head of the string. For the query example, this would be .before="select " . (Default: NULL) |
after | The text to put at the tail of the string. For the query example, .after=" from data_table" . (Default: NULL) |
between_cols | The text to insert between columns of text. See below for an example (Default is set to equal .between ) |
prune | If you don't want to use the entire text set, you can provide a function to indicate which elements should be pruned out. Some examples: |
prune_parameter | A void pointer to pass to your prune function. |
strings
table joined as per your specification. Allocated by the function, to be freed by you if desired.NULL
or has no text, I will print only the .before
and .after
parts with nothing in between. apop_opts.verbose >=3
, then print the pasted text to stderr. The sample snippet generates the SQL for a query using a list of column names (where the query begins with select
, ends with from datatab
, and has commas in between each element), re-processes the same list to produce the head of an HTML table, then produces the body of the table with the query result (pasting the tr
s and td
s into the right places).
double* apop_vector_percentiles | ( | gsl_vector * | data, |
char | rounding | ||
) |
Returns an array of size 101, where returned_vector
[95] gives the value of the 95th percentile, for example. Returned_vector
[100] is always the maximum value, and returned_vector
[0] is always the min (regardless of rounding rule).
data | a gsl_vector of data. (No default, must not be NULL .) |
rounding | This will either be 'u' , 'd' , or 'a' . Unless your data is exactly a multiple of 101, some percentiles will be ambiguous. If 'u' , then round up (use the next highest value); if 'd' (or anything else), round down to the next lowest value; if 'a' , take the mean of the two nearest points. If 'u' or 'a' , then you can say "5% or more of the sample is below \c returned_vector[5]"; if 'd' or 'a' , then you can say "5% or more of the sample is above returned_vector[5]". (Default = 'd' .) |
free()
the array returned by this function.