time-to-botec

Benchmark sampling in different programming languages
Log | Files | Refs | README

main.js (6055B)


      1 /**
      2 * @license Apache-2.0
      3 *
      4 * Copyright (c) 2018 The Stdlib Authors.
      5 *
      6 * Licensed under the Apache License, Version 2.0 (the "License");
      7 * you may not use this file except in compliance with the License.
      8 * You may obtain a copy of the License at
      9 *
     10 *    http://www.apache.org/licenses/LICENSE-2.0
     11 *
     12 * Unless required by applicable law or agreed to in writing, software
     13 * distributed under the License is distributed on an "AS IS" BASIS,
     14 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
     15 * See the License for the specific language governing permissions and
     16 * limitations under the License.
     17 */
     18 
     19 'use strict';
     20 
     21 // MODULES //
     22 
     23 var setReadOnly = require( '@stdlib/utils/define-read-only-property' );
     24 var isNumberArray = require( '@stdlib/assert/is-number-array' );
     25 var isNonNegativeInteger = require( '@stdlib/assert/is-nonnegative-integer' );
     26 var betaQuantile = require( './../../base/dists/beta/quantile' );
     27 var floor = require( '@stdlib/math/base/special/floor' );
     28 var ceil = require( '@stdlib/math/base/special/ceil' );
     29 var binomialCDF = require( './../../base/dists/binomial/cdf' );
     30 var binomialPMF = require( './../../base/dists/binomial/pmf' );
     31 var validate = require( './validate.js' );
     32 var print = require( './print.js' ); // eslint-disable-line stdlib/no-redeclare
     33 
     34 
     35 // VARIABLES //
     36 
     37 var RELATIVE_ERROR = 1+1e-07;
     38 
     39 
     40 // FUNCTIONS //
     41 
     42 /**
     43 * Calculates the lower endpoint of a confidence interval.
     44 *
     45 * @private
     46 * @param {NonNegativeInteger} x - number of successes
     47 * @param {NonNegativeInteger} n - total number of observations
     48 * @param {number} alpha - significance level
     49 * @returns {number} lower endpoint
     50 */
     51 function lower( x, n, alpha ) {
     52 	return ( x === 0 ) ? 0 : betaQuantile( alpha, x, n - x + 1 );
     53 }
     54 
     55 /**
     56 * Calculates the upper endpoint of a confidence interval.
     57 *
     58 * @private
     59 * @param {NonNegativeInteger} x - number of successes
     60 * @param {NonNegativeInteger} n - total number of observations
     61 * @param {number} alpha - significance level
     62 * @returns {number} upper endpoint
     63 */
     64 function upper( x, n, alpha ) {
     65 	return ( x === n ) ? 1 : betaQuantile( 1 - alpha, x + 1, n - x );
     66 }
     67 
     68 
     69 // MAIN //
     70 
     71 /**
     72 * Computes an exact test for the success probability in a Bernoulli experiment.
     73 *
     74 * @param {(NonNegativeInteger|Array)} x - number of successes or two-element array with successes and failures
     75 * @param {NonNegativeInteger} [n] - total number of observations
     76 * @param {Options} [options] - function options
     77 * @param {number} [options.alpha=0.05] - significance level
     78 * @param {string} [options.alternative='two-sided'] - alternative hypothesis (`two-sided`, `less` or `greater`)
     79 * @param {Probability} [options.p=0.5] - success probability under H0
     80 * @throws {TypeError} options has to be simple object
     81 * @throws {TypeError} must provide valid options
     82 * @throws {RangeError} alpha option has to be a number in the interval `[0,1]`
     83 * @throws {TypeError} alternative option has to be a string primitive
     84 * @throws {Error} alternative option must be `two-sided`, `less` or `greater`
     85 * @returns {Object} test results
     86 */
     87 function binomialTest() {
     88 	var alpha;
     89 	var opts;
     90 	var cint;
     91 	var pval;
     92 	var stat;
     93 	var alt;
     94 	var err;
     95 	var out;
     96 	var d;
     97 	var m;
     98 	var n;
     99 	var p;
    100 	var x;
    101 	var y;
    102 	var i;
    103 
    104 	opts = {};
    105 	if ( isNumberArray( arguments[ 0 ] ) ) {
    106 		x = arguments[ 0 ];
    107 		if ( x.length !== 2 ) {
    108 			throw new Error( 'invalid argument. If provided an array, it must have two elements. Value: `' + x + '`.' );
    109 		}
    110 		n = x[ 1 ] + x[ 0 ];
    111 		x = x[ 0 ];
    112 		if ( arguments[ 1 ] ) {
    113 			err = validate( opts, arguments[ 1 ] );
    114 		}
    115 	} else {
    116 		x = arguments[ 0 ];
    117 		n = arguments[ 1 ];
    118 		if ( !isNonNegativeInteger( x ) ) {
    119 			throw new TypeError( 'invalid argument. Must provide a nonnegative integer or a two-element array. Value: `' + x + '`.' );
    120 		}
    121 		if ( !isNonNegativeInteger( n ) ) {
    122 			throw new TypeError( 'invalid argument. Must provide a nonnegative integer. Value: `' + n + '`.' );
    123 		}
    124 		if ( x > n ) {
    125 			throw new TypeError( 'invalid arguments. `x` cannot be larger than `n`. `x:' + x + ', n:' + n + '`.' );
    126 		}
    127 		if ( arguments[ 2 ] ) {
    128 			err = validate( opts, arguments[ 2 ] );
    129 		}
    130 	}
    131 	if ( err ) {
    132 		throw err;
    133 	}
    134 
    135 	if ( opts.alpha === void 0 ) {
    136 		alpha = 0.05;
    137 	} else {
    138 		alpha = opts.alpha;
    139 	}
    140 	if ( alpha < 0.0 || alpha > 1.0 ) {
    141 		throw new RangeError( 'invalid argument. Option `alpha` must be a number in the range 0 to 1. Value: `' + alpha + '`.' );
    142 	}
    143 	if ( opts.p === void 0 ) {
    144 		p = 0.5;
    145 	} else {
    146 		p = opts.p;
    147 	}
    148 	if ( p < 0.0 || p > 1.0 ) {
    149 		throw new RangeError( 'invalid argument. Option `p` must be a probability. Value: `' + p + '`.' );
    150 	}
    151 
    152 	alt = opts.alternative || 'two-sided';
    153 	stat = x / n;
    154 	switch ( alt ) {
    155 	case 'less':
    156 		pval = binomialCDF( x, n, p );
    157 		cint = [ 0.0, upper( x, n, alpha ) ];
    158 		break;
    159 	case 'greater':
    160 		pval = 1.0 - binomialCDF( x - 1, n, p );
    161 		cint = [ lower( x, n, alpha ), 1.0 ];
    162 		break;
    163 	case 'two-sided':
    164 		d = binomialPMF( x, n, p );
    165 		m = n * p;
    166 		if ( x === m ) {
    167 			pval = 1;
    168 		} else if ( x < m ) {
    169 			y = 0;
    170 			for ( i = ceil( m ); i <= n; i++ ) {
    171 				if ( binomialPMF( i, n, p ) <= d * RELATIVE_ERROR ) {
    172 					y += 1;
    173 				}
    174 			}
    175 			pval = binomialCDF(x, n, p) + ( 1 - binomialCDF(n - y, n, p ) );
    176 		} else {
    177 			y = 0;
    178 			for ( i = 0; i <= floor( m ); i++ ) {
    179 				if ( binomialPMF( i, n, p ) <= d * RELATIVE_ERROR ) {
    180 					y += 1;
    181 				}
    182 			}
    183 			pval = binomialCDF( y-1, n, p ) + ( 1 - binomialCDF( x-1, n, p ) );
    184 		}
    185 		cint = [ lower( x, n, alpha/2.0 ), upper( x, n, alpha/2.0 ) ];
    186 		break;
    187 	default:
    188 		throw new Error( 'Invalid option. `alternative` must be either `two-sided`, `less` or `greater`. Value: `' + alt + '`' );
    189 	}
    190 
    191 	out = {};
    192 	setReadOnly( out, 'rejected', pval <= alpha );
    193 	setReadOnly( out, 'alpha', alpha );
    194 	setReadOnly( out, 'pValue', pval );
    195 	setReadOnly( out, 'statistic', stat );
    196 	setReadOnly( out, 'ci', cint );
    197 	setReadOnly( out, 'nullValue', p );
    198 	setReadOnly( out, 'alternative', alt );
    199 	setReadOnly( out, 'method', 'Exact binomial test' );
    200 	setReadOnly( out, 'print', print );
    201 	return out;
    202 }
    203 
    204 
    205 // EXPORTS //
    206 
    207 module.exports = binomialTest;