view release on metacpan or  search on metacpan

CaseResampling.xs  view on Meta::CPAN

double
mean(sample)
    AV* sample
  CODE:
    RETVAL = cs_mean_av(aTHX_ sample);
  OUTPUT: RETVAL


double
sample_standard_deviation(mean, sample)
    SV* mean
    AV* sample
  CODE:
    RETVAL =  cs_sum_deviation_squared_av(aTHX_ SvNV(mean), sample);
    RETVAL = pow( RETVAL / av_len(sample), 0.5 ); /* av_len() is N-1! */
  OUTPUT: RETVAL


double
population_standard_deviation(mean, sample)
    SV* mean
    AV* sample
  CODE:
    RETVAL =  cs_sum_deviation_squared_av(aTHX_ SvNV(mean), sample);
    RETVAL = pow( RETVAL / (av_len(sample)+1), 0.5 ); /* av_len() is N-1! */
  OUTPUT: RETVAL


double
select_kth(sample, kth)
    AV* sample
    I32 kth
  PREINIT:
    I32 nelem;
    double* csample;
  CODE:
    avToCAry(aTHX_ sample, &csample, &nelem);
    if (kth < 1 || kth > nelem) {
      croak("Can't select %ith smallest element from a list of %i elements", kth, nelem);
    }
    RETVAL = cs_select(csample, nelem, kth-1);
    Safefree(csample);
  OUTPUT: RETVAL


void
median_simple_confidence_limits(sample, confidence, ...)
    AV* sample
    double confidence
  PREINIT:
    /* "confidence" is 1-alpha */
    I32 runs, nelem, i_run;
    double *csample, *destsample, *medians;
    struct mt* rnd;
    double median   = 0.;
    double lower_ci = 0.;
    double upper_ci = 0.;
    double alpha;
  INIT:
    alpha = 1.-confidence;
  PPCODE:
    if (items == 2)
      runs = 1000;
    else if (items == 3)
      runs = SvUV(ST(2));
    else {
      croak("Usage: ($lower, $median, $upper) = median_confidence_limits(\\@sample, $confidence, [$nruns]);");
    }
    if (confidence <= 0. || confidence >= 1.) {
      croak("Confidence level has to be in (0, 1)");
    }
    rnd = get_rnd(aTHX);
    avToCAry(aTHX_ sample, &csample, &nelem);
    if (nelem != 0) {
      median = cs_median(csample, nelem);
      Newx(medians, runs, double);
      Newx(destsample, nelem, double);
      for (i_run = 0; i_run < runs; ++i_run) {
        do_resample(csample, nelem, rnd, destsample);
        medians[i_run] = cs_median(destsample, nelem);
      }
      Safefree(destsample);
      /* lower = t - (t*_((R+1)*(1-alpha)) - t)
       * upper = t - (t*_((R+1)*alpha) - t)
       */
      lower_ci = 2.*median - cs_select( medians, runs, (I32)((runs+1.)*(1.-alpha)) );
      upper_ci = 2.*median - cs_select( medians, runs, (I32)((runs+1.)*alpha) );
      Safefree(medians);
    }
    Safefree(csample);
    EXTEND(SP, 3);
    mPUSHn(lower_ci);
    mPUSHn(median);
    mPUSHn(upper_ci);


void
simple_confidence_limits_from_samples(statistic, statistics, confidence)
    double statistic
    AV* statistics
    double confidence
  PREINIT:
    /* "confidence" is 1-alpha */
    I32 nelem;
    double *cstatistics;
    double lower_ci = 0.;
    double upper_ci = 0.;
    double alpha;
  INIT:
    alpha = 1.-confidence;
  PPCODE:
    if (confidence <= 0. || confidence >= 1.) {
      croak("Confidence level has to be in (0, 1)");
    }
    avToCAry(aTHX_ statistics, &cstatistics, &nelem);
    if (nelem != 0) {
      /* lower = t - (t*_((R+1)*(1-alpha)) - t)
       * upper = t - (t*_((R+1)*alpha) - t)
       */
      lower_ci = 2.*statistic - cs_select( cstatistics, nelem, (I32)((nelem+1.)*(1.-alpha)) );
      upper_ci = 2.*statistic - cs_select( cstatistics, nelem, (I32)((nelem+1.)*alpha) );
    }
    Safefree(cstatistics);
    EXTEND(SP, 3);
    mPUSHn(lower_ci);
    mPUSHn(statistic);
    mPUSHn(upper_ci);


double
approx_erf(x)
    double x
  CODE:
    RETVAL = cs_approx_erf(x);
  OUTPUT: RETVAL


double
approx_erf_inv(x)
    double x
  CODE:
    if (x <= 0. || x >= 1.)
      croak("The inverse error function is defined in (0,1). %f is outside that range", x);
    RETVAL = cs_approx_erf_inv(x);
  OUTPUT: RETVAL


double
alpha_to_nsigma(x)
    double x
  CODE:
    RETVAL = cs_alpha_to_nsigma(x);
  OUTPUT: RETVAL


double
nsigma_to_alpha(x)
    double x
  CODE:
    RETVAL = cs_nsigma_to_alpha(x);
  OUTPUT: RETVAL