time-to-botec

kde2d.js (5295B)

---

 /**
 * @license Apache-2.0
 *
 * Copyright (c) 2018 The Stdlib Authors.
 *
 * Licensed under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License.
 * You may obtain a copy of the License at
 *
 *    http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
 
 'use strict';
 
 // MODULES //
 
 var ndarray = require( '@stdlib/ndarray/array' );
 var linspace = require( '@stdlib/array/linspace' );
 var setReadOnly = require( '@stdlib/utils/define-read-only-property' );
 var isNumericArray = require( '@stdlib/assert/is-numeric-array' );
 var isMatrixLike = require( '@stdlib/assert/is-matrix-like' );
 var pickBandwidth = require( './pick_bandwidth.js' );
 var validate = require( './validate.js' );
 var ndarrayLike = require( './ndarray_like.js' );
 var min = require( './min.js' );
 var max = require( './max.js' );
 var gaussian = require( './gaussian.js' );
 
 
 // MAIN //
 
 /**
 * Computes two-dimensional kernel density estimates.
 *
 * @param {NumericArray} x - array of x values
 * @param {NumericArray} y - array of y values
 * @param {Options} [options] - function options
 * @param {NumericArray} [options.h] - array of length two containing the bandwidth values for x and y
 * @param {number} [options.n=25] - number of partitions on the x- and y-axes
 * @param {number} [options.xMin] - lower limit of x
 * @param {number} [options.xMax] - upper limit of x
 * @param {number} [options.yMin] - lower limit of y
 * @param {number} [options.yMax] - upper limit of y
 * @param {(string|Function)} [options.kernel='gaussian'] - a string or function to specifying the used kernel function
 * @throws {TypeError} first argument must be an array or matrix-like
 * @throws {TypeError} second argument must be an array
 * @throws {Error} first and second arguments must be of the same length
 * @throws {RangeError} `xMin` must be smaller than `xMax`
 * @throws {RangeError} `yMin` must be smaller than `yMax`
 * @throws {TypeError} options argument must be an object
 * @throws {TypeError} must provide valid options
 * @returns {Object} object containing the density estimates (`z`) along grid points (`x` and `y` values)
 *
 * @example
 * var x = [ 0.6333, 0.8643, 1.0952, 1.3262, 1.5571, 1.7881, 2.019, 2.25, 2.481, 2.7119 ];
 * var y = [ -0.0468, 0.8012, 1.6492, 2.4973, 3.3454, 4.1934, 5.0415, 5.8896, 6.7376, 7.5857 ];
 * var out = kde2d( x, y );
 */
 function kde2d() {
 	var kernelFunction;
 	var maxArgs;
 	var zScoreX;
 	var zScoreY;
 	var gridX;
 	var gridY;
 	var xMin;
 	var xMax;
 	var yMin;
 	var yMax;
 	var xVal; // For gridspace loop
 	var yVal; // For gridspace loop
 	var subX;
 	var subY;
 	var opts;
 	var arr;
 	var err;
 	var ans;
 	var out;
 	var gx;
 	var gy;
 	var hX;
 	var hY;
 	var ix;
 	var iy;
 	var x;
 	var y;
 	var i;
 	var n;
 	var z;
 
 	opts = {};
 
 	if ( isMatrixLike( arguments[0] ) ) {
 		// Case of ndarray, opts
 		arr = arguments[ 0 ];
 		maxArgs = 1;
 	} else {
 		x = arguments[ 0 ];
 		y = arguments[ 1 ];
 		if ( !isNumericArray( x ) ) {
 			throw new TypeError( 'invalid argument. First argument `x` must be a numeric array. Value: `' + x + '`.' );
 		}
 		if ( !isNumericArray( y ) ) {
 			throw new TypeError( 'invalid argument. Second argument `y` must be a numeric array. Value: `' + y + '`.' );
 		}
 		if ( x.length !== y.length ) {
 			throw new Error( 'invalid arguments. Arguments `x` and `y` must be arrays of the same length' );
 		}
 		arr = ndarrayLike( x, y );
 		maxArgs = 2;
 	}
 
 	if ( arguments.length > maxArgs ) {
 		err = validate( opts, arguments[ maxArgs ] );
 		if ( err ) {
 			throw err;
 		}
 	}
 
 	if ( opts.h ) {
 		hX = opts.h[0];
 		hY = opts.h[1];
 	} else {
 		hX = pickBandwidth(arr, 0);
 		hY = pickBandwidth(arr, 1);
 	}
 
 	n = opts.n || 25;
 	xMin = opts.xMin || min( arr, 0, arr.shape[0] );
 	xMax = opts.xMax || max( arr, 0, arr.shape[0] );
 	yMin = opts.yMin || min( arr, 1, arr.shape[0] );
 	yMax = opts.yMax || max( arr, 1, arr.shape[0] );
 
 	if ( xMin >= xMax ) {
 		throw new RangeError( '`x` min must be strictly less than max' );
 	}
 	if ( yMin >= yMax ) {
 		throw new RangeError( '`y` min must be strictly less than max' );
 	}
 
 	kernelFunction = opts.kernel || gaussian;
 
 	// Create the `ndarray` to hold the density values:
 	z = ndarray({
 		'shape': [n, n]
 	} );
 
 	// Make the grid:
 	gridX = linspace(xMin, xMax, n);
 	gridY = linspace(yMin, yMax, n);
 
 	// Loop through x and y indices:
 	for ( ix = 0; ix < gridX.length; ix++ ) {
 		gx = gridX[ ix ];
 		for ( iy = 0; iy < gridY.length; iy++ ) {
 			gy = gridY[ iy ];
 			ans = 0.0;
 			for ( i = 0; i < arr.shape[ 0 ]; i++ ) {
 				xVal = arr.get( i, 0 );
 				yVal = arr.get( i, 1 );
 
 				zScoreX = ( (xVal - gx) / hX );
 				zScoreY = ( (yVal - gy) / hY );
 
 				subX = ( 1.0 / hX ) * kernelFunction( zScoreX );
 				subY = ( 1.0 / hY ) * kernelFunction( zScoreY );
 				ans += ( subX * subY );
 			}
 			z.set( ix, iy, ans / arr.shape[0] );
 		}
 	}
 	out = {};
 	setReadOnly( out, 'x', gridX );
 	setReadOnly( out, 'y', gridY );
 	setReadOnly( out, 'z', z );
 	return out;
 }
 
 
 // EXPORTS //
 
 module.exports = kde2d;