time-to-botec

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

main.js (5635B)


      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 isCollection = require( '@stdlib/assert/is-collection' );
     24 var isPlainObject = require( '@stdlib/assert/is-plain-object' );
     25 var setReadOnly = require( '@stdlib/utils/define-read-only-property' );
     26 var countBy = require( '@stdlib/utils/count-by' );
     27 var objectKeys = require( '@stdlib/utils/keys' );
     28 var rank = require( './../../ranks' );
     29 var pow = require( '@stdlib/math/base/special/pow' );
     30 var chisqCDF = require( './../../base/dists/chisquare/cdf' );
     31 var identity = require( '@stdlib/utils/identity-function' );
     32 var incrspace = require( '@stdlib/array/incrspace' );
     33 var validate = require( './validate.js' );
     34 var print = require( './print.js' ); // eslint-disable-line stdlib/no-redeclare
     35 
     36 
     37 // MAIN //
     38 
     39 /**
     40 * Computes the Kruskal-Wallis test for equality of medians.
     41 *
     42 * @param {...NumberArray} arguments - either two or more number arrays or a single numeric array if an array of group indicators is supplied as an option
     43 * @param {Options} [options] - function options
     44 * @param {number} [options.alpha=0.05] - significance level
     45 * @param {Array} [options.groups] - array of group indicators
     46 * @throws {Error} must provide at least two array-like arguments if `groups` is not set
     47 * @throws {TypeError} must provide array-like arguments
     48 * @throws {TypeError} options has to be simple object
     49 * @throws {TypeError} must provide valid options
     50 * @throws {RangeError} alpha option has to be a number in the interval `[0,1]`
     51 * @returns {Object} test results
     52 *
     53 * @example
     54 * // Data from Hollander & Wolfe (1973), p. 116:
     55 * var x = [ 2.9, 3.0, 2.5, 2.6, 3.2 ];
     56 * var y = [ 3.8, 2.7, 4.0, 2.4 ];
     57 * var z = [ 2.8, 3.4, 3.7, 2.2, 2.0 ];
     58 *
     59 * var out = kruskal( x, y, z );
     60 * // returns {...}
     61 */
     62 function kruskal() {
     63 	var groupsIndicators;
     64 	var groupRankSums;
     65 	var tieSumTerm;
     66 	var ngroups;
     67 	var options;
     68 	var levels;
     69 	var alpha;
     70 	var param;
     71 	var ranks;
     72 	var vals;
     73 	var opts;
     74 	var pval;
     75 	var stat;
     76 	var ties;
     77 	var arg;
     78 	var err;
     79 	var key;
     80 	var out;
     81 	var i;
     82 	var j;
     83 	var n;
     84 	var N;
     85 	var x;
     86 	var v;
     87 
     88 	ngroups = arguments.length;
     89 	opts = {};
     90 	if ( isPlainObject( arguments[ ngroups - 1 ] ) ) {
     91 		options = arguments[ ngroups - 1 ];
     92 		ngroups -= 1;
     93 		err = validate( opts, options );
     94 		if ( err ) {
     95 			throw err;
     96 		}
     97 	}
     98 	groupRankSums = {};
     99 	n = {};
    100 	if ( opts.groups ) {
    101 		x = arguments[ 0 ];
    102 		if ( x.length !== opts.groups.length ) {
    103 			throw new RangeError( 'invalid arguments. First argument and `opts.groups` must be arrays of the same length.' );
    104 		}
    105 		n = countBy( opts.groups, identity );
    106 		levels = objectKeys( n );
    107 		ngroups = levels.length;
    108 		for ( i = 0; i < ngroups; i++ ) {
    109 			key = levels[ i ];
    110 			groupRankSums[ key ] = 0;
    111 		}
    112 		if ( ngroups < 2 ) {
    113 			throw new Error( 'invalid number of groups. `groups` array must contain at least two unique elements. Value: `' + levels + '`.' );
    114 		}
    115 		groupsIndicators = opts.groups;
    116 	} else {
    117 		x = [];
    118 		groupsIndicators = [];
    119 		if ( ngroups < 2 ) {
    120 			throw new Error( 'invalid number of input arguments. Must provide at least two array-like arguments. Value: `' + arg + '`.' );
    121 		}
    122 		for ( i = 0; i < ngroups; i++ ) {
    123 			arg = arguments[ i ];
    124 			if ( !isCollection( arg ) ) {
    125 				throw new TypeError( 'invalid argument. Must provide array-like arguments. Value: `' + arg + '`.' );
    126 			}
    127 			if ( arg.length === 0 ) {
    128 				throw new Error( 'invalid argument. Supplied arrays cannot be empty. Value: `' + arg + '`.' );
    129 			} else {
    130 				n[ i ] = arg.length;
    131 			}
    132 			groupRankSums[ i ] = 0;
    133 			for ( j = 0; j < n[ i ]; j++ ) {
    134 				groupsIndicators.push( i );
    135 				x.push( arg[ j ] );
    136 			}
    137 		}
    138 		levels = incrspace( 0, ngroups, 1 );
    139 	}
    140 	if ( opts.alpha === void 0 ) {
    141 		alpha = 0.05;
    142 	} else {
    143 		alpha = opts.alpha;
    144 	}
    145 	if ( alpha < 0.0 || alpha > 1.0 ) {
    146 		throw new RangeError( 'invalid option. `alpha` must be a number in the range 0 to 1. Value: `' + alpha + '`.' );
    147 	}
    148 
    149 	N = x.length;
    150 	ranks = rank( x );
    151 
    152 	// Calculate # ties for each value & rank sums per group:
    153 	ties = {};
    154 	for ( i = 0; i < N; i++ ) {
    155 		groupRankSums[ groupsIndicators[ i ] ] += ranks[ i ];
    156 		if ( x[ i ] in ties ) {
    157 			ties[ x[ i ] ] += 1;
    158 		} else {
    159 			ties[ x[ i ] ] = 1;
    160 		}
    161 	}
    162 
    163 	// Calculate test statistic using short-cut formula:
    164 	stat = 0.0;
    165 	for ( i = 0; i < ngroups; i++ ) {
    166 		key = levels[ i ];
    167 		stat += pow( groupRankSums[ key ], 2.0 ) / n[ key ];
    168 	}
    169 	stat = ( ( 12.0 / ( N * (N+1) ) ) * stat ) - ( 3.0 * (N+1) );
    170 
    171 	// Correction for ties:
    172 	tieSumTerm = 0;
    173 	vals = objectKeys( ties );
    174 	for ( i = 0; i < vals.length; i++ ) {
    175 		v = ties[ vals[ i ] ];
    176 		tieSumTerm += pow( v, 3.0 ) - v;
    177 	}
    178 
    179 	stat /= 1.0 - ( ( tieSumTerm ) / ( pow( N, 3 ) - N ) );
    180 	param = ngroups - 1;
    181 	pval = 1.0 - chisqCDF( stat, param );
    182 
    183 	out = {};
    184 	setReadOnly( out, 'rejected', pval <= alpha );
    185 	setReadOnly( out, 'alpha', alpha );
    186 	setReadOnly( out, 'df', param );
    187 	setReadOnly( out, 'pValue', pval );
    188 	setReadOnly( out, 'statistic', stat );
    189 	setReadOnly( out, 'method', 'Kruskal-Wallis Test' );
    190 	setReadOnly( out, 'print', print );
    191 	return out;
    192 }
    193 
    194 
    195 // EXPORTS //
    196 
    197 module.exports = kruskal;