Giu/analistica
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

ex-1: reduce sample size and various fixes

Browse Source
This commit is contained in:
Michele Guerini Rocco 2020-04-11 19:29:38 +02:00
parent 72fa3ad2e7
commit 343fe3aac3
3 changed files with 7 additions and 7 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
ex-1/bootstrap.c
View File

@ -43,8 +43,8 @@ size_t mean_index(gsl_vector *v, double x) {
/* Computes the half-sample mode (also called the Robertson-Cryer /* Computes the half-sample mode (also called the Robertson-Cryer
* mode estimator) of the sample `x` containing `n` observations. * mode estimator) of the sample `x` containing `n` observations.
* *
* It is based on repeatedly finding the modal interval (interval * It is based on repeatedly finding the half modal interval (smallest
* containing the most observations) of half of the sample. * interval containing half of the observations) of the sample.
* This implementation is based on the `hsm()` function from the * This implementation is based on the `hsm()` function from the
* modeest[1] R package. * modeest[1] R package.
* *
@ -129,7 +129,7 @@ double gauss_kde(double x, void * params) {
* by the sample variance times a factor which * by the sample variance times a factor which
* depends on the number of points and dimension. * depends on the number of points and dimension.
*/ */
double bw = p.var * pow((double)p.n*3.0/4, -2.0/5); double bw = 0.4 * p.var * pow((double)p.n*3.0/4, -2.0/5);
double sum = 0; double sum = 0;
for (size_t i = 0; i < p.n; i++) for (size_t i = 0; i < p.n; i++)

6
ex-1/main.c
View File

@ -20,7 +20,7 @@ int main(int argc, char** argv) {
gsl_rng *r = gsl_rng_alloc(gsl_rng_default); gsl_rng *r = gsl_rng_alloc(gsl_rng_default);
// prepare histogram // prepare histogram
size_t samples = 100000; size_t samples = 50000;
double* sample = calloc(samples, sizeof(double)); double* sample = calloc(samples, sizeof(double));
double min = -10; double min = -10;
double max = 10; double max = 10;
@ -40,7 +40,7 @@ int main(int argc, char** argv) {
} }
fprintf(stderr, "done\n"); fprintf(stderr, "done\n");
// sort the sample: needed for HSM and KS tests // sort the sample: needed for HSM and ks tests
qsort(sample, samples, sizeof(double), &cmp_double); qsort(sample, samples, sizeof(double), &cmp_double);
@ -88,7 +88,7 @@ int main(int argc, char** argv) {
// print the results // print the results
fprintf(stderr, "\n## Results\n"); fprintf(stderr, "\n## Results\n");
fprintf(stderr, "expected mode: %.7f\n", mode_e); fprintf(stderr, "expected mode: %.8f\n", mode_e);
fprintf(stderr, "observed mode: %.4f±%.4f\n", mode_o.n, mode_o.s); fprintf(stderr, "observed mode: %.4f±%.4f\n", mode_o.n, mode_o.s);
// t-test // t-test

2
ex-1/tests.c
View File

@ -108,7 +108,7 @@ double numeric_mode(double min, double max,
x = gsl_min_fminimizer_x_minimum(s); x = gsl_min_fminimizer_x_minimum(s);
min = gsl_min_fminimizer_x_lower(s); min = gsl_min_fminimizer_x_lower(s);
max = gsl_min_fminimizer_x_upper(s); max = gsl_min_fminimizer_x_upper(s);
status = gsl_min_test_interval(min, max, prec, prec); status = gsl_min_test_interval(min, max, 0, prec);
} while (status == GSL_CONTINUE && iter < max_iter); } while (status == GSL_CONTINUE && iter < max_iter);