| /* |
| * The MIT License (MIT) |
| * |
| * Copyright (c) 2015 ml.js |
| * |
| * Permission is hereby granted, free of charge, to any person obtaining a copy |
| * of this software and associated documentation files (the "Software"), to deal |
| * in the Software without restriction, including without limitation the rights |
| * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell |
| * copies of the Software, and to permit persons to whom the Software is |
| * furnished to do so, subject to the following conditions: |
| * |
| * The above copyright notice and this permission notice shall be included in all |
| * copies or substantial portions of the Software. |
| * |
| * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR |
| * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, |
| * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE |
| * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER |
| * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, |
| * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE |
| * SOFTWARE. |
| */ |
| 'use strict'; |
| |
| // ml-stat array.js |
| const MLStatArray = {}; |
| { |
| function compareNumbers(a, b) { |
| return a - b; |
| } |
| |
| /** |
| * Computes the sum of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.sum = function sum(values) { |
| var sum = 0; |
| for (var i = 0; i < values.length; i++) { |
| sum += values[i]; |
| } |
| return sum; |
| }; |
| |
| /** |
| * Computes the maximum of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.max = function max(values) { |
| var max = values[0]; |
| var l = values.length; |
| for (var i = 1; i < l; i++) { |
| if (values[i] > max) max = values[i]; |
| } |
| return max; |
| }; |
| |
| /** |
| * Computes the minimum of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.min = function min(values) { |
| var min = values[0]; |
| var l = values.length; |
| for (var i = 1; i < l; i++) { |
| if (values[i] < min) min = values[i]; |
| } |
| return min; |
| }; |
| |
| /** |
| * Computes the min and max of the given values |
| * @param {Array} values |
| * @returns {{min: number, max: number}} |
| */ |
| MLStatArray.minMax = function minMax(values) { |
| var min = values[0]; |
| var max = values[0]; |
| var l = values.length; |
| for (var i = 1; i < l; i++) { |
| if (values[i] < min) min = values[i]; |
| if (values[i] > max) max = values[i]; |
| } |
| return { |
| min: min, |
| max: max |
| }; |
| }; |
| |
| /** |
| * Computes the arithmetic mean of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.arithmeticMean = function arithmeticMean(values) { |
| var sum = 0; |
| var l = values.length; |
| for (var i = 0; i < l; i++) { |
| sum += values[i]; |
| } |
| return sum / l; |
| }; |
| |
| /** |
| * {@link arithmeticMean} |
| */ |
| MLStatArray.mean = MLStatArray.arithmeticMean; |
| |
| /** |
| * Computes the geometric mean of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.geometricMean = function geometricMean(values) { |
| var mul = 1; |
| var l = values.length; |
| for (var i = 0; i < l; i++) { |
| mul *= values[i]; |
| } |
| return Math.pow(mul, 1 / l); |
| }; |
| |
| /** |
| * Computes the mean of the log of the given values |
| * If the return value is exponentiated, it gives the same result as the |
| * geometric mean. |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.logMean = function logMean(values) { |
| var lnsum = 0; |
| var l = values.length; |
| for (var i = 0; i < l; i++) { |
| lnsum += Math.log(values[i]); |
| } |
| return lnsum / l; |
| }; |
| |
| /** |
| * Computes the weighted grand mean for a list of means and sample sizes |
| * @param {Array} means - Mean values for each set of samples |
| * @param {Array} samples - Number of original values for each set of samples |
| * @returns {number} |
| */ |
| MLStatArray.grandMean = function grandMean(means, samples) { |
| var sum = 0; |
| var n = 0; |
| var l = means.length; |
| for (var i = 0; i < l; i++) { |
| sum += samples[i] * means[i]; |
| n += samples[i]; |
| } |
| return sum / n; |
| }; |
| |
| /** |
| * Computes the truncated mean of the given values using a given percentage |
| * @param {Array} values |
| * @param {number} percent - The percentage of values to keep (range: [0,1]) |
| * @param {boolean} [alreadySorted=false] |
| * @returns {number} |
| */ |
| MLStatArray.truncatedMean = function truncatedMean(values, percent, alreadySorted) { |
| if (alreadySorted === undefined) alreadySorted = false; |
| if (!alreadySorted) { |
| values = [].concat(values).sort(compareNumbers); |
| } |
| var l = values.length; |
| var k = Math.floor(l * percent); |
| var sum = 0; |
| for (var i = k; i < (l - k); i++) { |
| sum += values[i]; |
| } |
| return sum / (l - 2 * k); |
| }; |
| |
| /** |
| * Computes the harmonic mean of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.harmonicMean = function harmonicMean(values) { |
| var sum = 0; |
| var l = values.length; |
| for (var i = 0; i < l; i++) { |
| if (values[i] === 0) { |
| throw new RangeError('value at index ' + i + 'is zero'); |
| } |
| sum += 1 / values[i]; |
| } |
| return l / sum; |
| }; |
| |
| /** |
| * Computes the contraharmonic mean of the given values |
| * @param {Array} values |
| * @returns {number} |
| */ |
| MLStatArray.contraHarmonicMean = function contraHarmonicMean(values) { |
| var r1 = 0; |
| var r2 = 0; |
| var l = values.length; |
| for (var i = 0; i < l; i++) { |
| r1 += values[i] * values[i]; |
| r2 += values[i]; |
| } |
| if (r2 < 0) { |
| throw new RangeError('sum of values is negative'); |
| } |
| return r1 / r2; |
| }; |
| |
| /** |
| * Computes the median of the given values |
| * @param {Array} values |
| * @param {boolean} [alreadySorted=false] |
| * @returns {number} |
| */ |
| MLStatArray.median = function median(values, alreadySorted) { |
| if (alreadySorted === undefined) alreadySorted = false; |
| if (!alreadySorted) { |
| values = [].concat(values).sort(compareNumbers); |
| } |
| var l = values.length; |
| var half = Math.floor(l / 2); |
| if (l % 2 === 0) { |
| return (values[half - 1] + values[half]) * 0.5; |
| } else { |
| return values[half]; |
| } |
| }; |
| |
| /** |
| * Computes the variance of the given values |
| * @param {Array} values |
| * @param {boolean} [unbiased=true] - if true, divide by (n-1); if false, divide by n. |
| * @returns {number} |
| */ |
| MLStatArray.variance = function variance(values, unbiased) { |
| if (unbiased === undefined) unbiased = true; |
| var theMean = MLStatArray.mean(values); |
| var theVariance = 0; |
| var l = values.length; |
| |
| for (var i = 0; i < l; i++) { |
| var x = values[i] - theMean; |
| theVariance += x * x; |
| } |
| |
| if (unbiased) { |
| return theVariance / (l - 1); |
| } else { |
| return theVariance / l; |
| } |
| }; |
| |
| /** |
| * Computes the standard deviation of the given values |
| * @param {Array} values |
| * @param {boolean} [unbiased=true] - if true, divide by (n-1); if false, divide by n. |
| * @returns {number} |
| */ |
| MLStatArray.standardDeviation = function standardDeviation(values, unbiased) { |
| return Math.sqrt(MLStatArray.variance(values, unbiased)); |
| }; |
| |
| MLStatArray.standardError = function standardError(values) { |
| return MLStatArray.standardDeviation(values) / Math.sqrt(values.length); |
| }; |
| |
| /** |
| * IEEE Transactions on biomedical engineering, vol. 52, no. 1, january 2005, p. 76- |
| * Calculate the standard deviation via the Median of the absolute deviation |
| * The formula for the standard deviation only holds for Gaussian random variables. |
| * @returns {{mean: number, stdev: number}} |
| */ |
| MLStatArray.robustMeanAndStdev = function robustMeanAndStdev(y) { |
| var mean = 0, stdev = 0; |
| var length = y.length, i = 0; |
| for (i = 0; i < length; i++) { |
| mean += y[i]; |
| } |
| mean /= length; |
| var averageDeviations = new Array(length); |
| for (i = 0; i < length; i++) |
| averageDeviations[i] = Math.abs(y[i] - mean); |
| averageDeviations.sort(compareNumbers); |
| if (length % 2 === 1) { |
| stdev = averageDeviations[(length - 1) / 2] / 0.6745; |
| } else { |
| stdev = 0.5 * (averageDeviations[length / 2] + averageDeviations[length / 2 - 1]) / 0.6745; |
| } |
| |
| return { |
| mean: mean, |
| stdev: stdev |
| }; |
| }; |
| |
| MLStatArray.quartiles = function quartiles(values, alreadySorted) { |
| if (typeof (alreadySorted) === 'undefined') alreadySorted = false; |
| if (!alreadySorted) { |
| values = [].concat(values).sort(compareNumbers); |
| } |
| |
| var quart = values.length / 4; |
| var q1 = values[Math.ceil(quart) - 1]; |
| var q2 = MLStatArray.median(values, true); |
| var q3 = values[Math.ceil(quart * 3) - 1]; |
| |
| return {q1: q1, q2: q2, q3: q3}; |
| }; |
| |
| MLStatArray.pooledStandardDeviation = function pooledStandardDeviation(samples, unbiased) { |
| return Math.sqrt(MLStatArray.pooledVariance(samples, unbiased)); |
| }; |
| |
| MLStatArray.pooledVariance = function pooledVariance(samples, unbiased) { |
| if (typeof (unbiased) === 'undefined') unbiased = true; |
| var sum = 0; |
| var length = 0, l = samples.length; |
| for (var i = 0; i < l; i++) { |
| var values = samples[i]; |
| var vari = MLStatArray.variance(values); |
| |
| sum += (values.length - 1) * vari; |
| |
| if (unbiased) |
| length += values.length - 1; |
| else |
| length += values.length; |
| } |
| return sum / length; |
| }; |
| |
| MLStatArray.mode = function mode(values) { |
| var l = values.length, |
| itemCount = new Array(l), |
| i; |
| for (i = 0; i < l; i++) { |
| itemCount[i] = 0; |
| } |
| var itemArray = new Array(l); |
| var count = 0; |
| |
| for (i = 0; i < l; i++) { |
| var index = itemArray.indexOf(values[i]); |
| if (index >= 0) |
| itemCount[index]++; |
| else { |
| itemArray[count] = values[i]; |
| itemCount[count] = 1; |
| count++; |
| } |
| } |
| |
| var maxValue = 0, maxIndex = 0; |
| for (i = 0; i < count; i++) { |
| if (itemCount[i] > maxValue) { |
| maxValue = itemCount[i]; |
| maxIndex = i; |
| } |
| } |
| |
| return itemArray[maxIndex]; |
| }; |
| |
| MLStatArray.covariance = function covariance(vector1, vector2, unbiased) { |
| if (typeof (unbiased) === 'undefined') unbiased = true; |
| var mean1 = MLStatArray.mean(vector1); |
| var mean2 = MLStatArray.mean(vector2); |
| |
| if (vector1.length !== vector2.length) |
| throw 'Vectors do not have the same dimensions'; |
| |
| var cov = 0, l = vector1.length; |
| for (var i = 0; i < l; i++) { |
| var x = vector1[i] - mean1; |
| var y = vector2[i] - mean2; |
| cov += x * y; |
| } |
| |
| if (unbiased) |
| return cov / (l - 1); |
| else |
| return cov / l; |
| }; |
| |
| MLStatArray.skewness = function skewness(values, unbiased) { |
| if (typeof (unbiased) === 'undefined') unbiased = true; |
| var theMean = MLStatArray.mean(values); |
| |
| var s2 = 0, s3 = 0, l = values.length; |
| for (var i = 0; i < l; i++) { |
| var dev = values[i] - theMean; |
| s2 += dev * dev; |
| s3 += dev * dev * dev; |
| } |
| var m2 = s2 / l; |
| var m3 = s3 / l; |
| |
| var g = m3 / (Math.pow(m2, 3 / 2.0)); |
| if (unbiased) { |
| var a = Math.sqrt(l * (l - 1)); |
| var b = l - 2; |
| return (a / b) * g; |
| } else { |
| return g; |
| } |
| }; |
| |
| MLStatArray.kurtosis = function kurtosis(values, unbiased) { |
| if (typeof (unbiased) === 'undefined') unbiased = true; |
| var theMean = MLStatArray.mean(values); |
| var n = values.length, s2 = 0, s4 = 0; |
| |
| for (var i = 0; i < n; i++) { |
| var dev = values[i] - theMean; |
| s2 += dev * dev; |
| s4 += dev * dev * dev * dev; |
| } |
| var m2 = s2 / n; |
| var m4 = s4 / n; |
| |
| if (unbiased) { |
| var v = s2 / (n - 1); |
| var a = (n * (n + 1)) / ((n - 1) * (n - 2) * (n - 3)); |
| var b = s4 / (v * v); |
| var c = ((n - 1) * (n - 1)) / ((n - 2) * (n - 3)); |
| |
| return a * b - 3 * c; |
| } else { |
| return m4 / (m2 * m2) - 3; |
| } |
| }; |
| |
| MLStatArray.entropy = function entropy(values, eps) { |
| if (typeof (eps) === 'undefined') eps = 0; |
| var sum = 0, l = values.length; |
| for (var i = 0; i < l; i++) |
| sum += values[i] * Math.log(values[i] + eps); |
| return -sum; |
| }; |
| |
| MLStatArray.weightedMean = function weightedMean(values, weights) { |
| var sum = 0, l = values.length; |
| for (var i = 0; i < l; i++) |
| sum += values[i] * weights[i]; |
| return sum; |
| }; |
| |
| MLStatArray.weightedStandardDeviation = function weightedStandardDeviation(values, weights) { |
| return Math.sqrt(MLStatArray.weightedVariance(values, weights)); |
| }; |
| |
| MLStatArray.weightedVariance = function weightedVariance(values, weights) { |
| var theMean = MLStatArray.weightedMean(values, weights); |
| var vari = 0, l = values.length; |
| var a = 0, b = 0; |
| |
| for (var i = 0; i < l; i++) { |
| var z = values[i] - theMean; |
| var w = weights[i]; |
| |
| vari += w * (z * z); |
| b += w; |
| a += w * w; |
| } |
| |
| return vari * (b / (b * b - a)); |
| }; |
| |
| MLStatArray.center = function center(values, inPlace) { |
| if (typeof (inPlace) === 'undefined') inPlace = false; |
| |
| var result = values; |
| if (!inPlace) |
| result = [].concat(values); |
| |
| var theMean = MLStatArray.mean(result), l = result.length; |
| for (var i = 0; i < l; i++) |
| result[i] -= theMean; |
| }; |
| |
| MLStatArray.standardize = function standardize(values, standardDev, inPlace) { |
| if (typeof (standardDev) === 'undefined') standardDev = MLStatArray.standardDeviation(values); |
| if (typeof (inPlace) === 'undefined') inPlace = false; |
| var l = values.length; |
| var result = inPlace ? values : new Array(l); |
| for (var i = 0; i < l; i++) |
| result[i] = values[i] / standardDev; |
| return result; |
| }; |
| |
| MLStatArray.cumulativeSum = function cumulativeSum(array) { |
| var l = array.length; |
| var result = new Array(l); |
| result[0] = array[0]; |
| for (var i = 1; i < l; i++) |
| result[i] = result[i - 1] + array[i]; |
| return result; |
| }; |
| } |
| |
| |
| // ml-stat matrix.js |
| const MLStatMatrix = {}; |
| { |
| let arrayStat = MLStatArray; |
| |
| function compareNumbers(a, b) { |
| return a - b; |
| } |
| |
| MLStatMatrix.max = function max(matrix) { |
| var max = -Infinity; |
| for (var i = 0; i < matrix.length; i++) { |
| for (var j = 0; j < matrix[i].length; j++) { |
| if (matrix[i][j] > max) max = matrix[i][j]; |
| } |
| } |
| return max; |
| }; |
| |
| MLStatMatrix.min = function min(matrix) { |
| var min = Infinity; |
| for (var i = 0; i < matrix.length; i++) { |
| for (var j = 0; j < matrix[i].length; j++) { |
| if (matrix[i][j] < min) min = matrix[i][j]; |
| } |
| } |
| return min; |
| }; |
| |
| MLStatMatrix.minMax = function minMax(matrix) { |
| var min = Infinity; |
| var max = -Infinity; |
| for (var i = 0; i < matrix.length; i++) { |
| for (var j = 0; j < matrix[i].length; j++) { |
| if (matrix[i][j] < min) min = matrix[i][j]; |
| if (matrix[i][j] > max) max = matrix[i][j]; |
| } |
| } |
| return { |
| min:min, |
| max:max |
| }; |
| }; |
| |
| MLStatMatrix.entropy = function entropy(matrix, eps) { |
| if (typeof (eps) === 'undefined') { |
| eps = 0; |
| } |
| var sum = 0, |
| l1 = matrix.length, |
| l2 = matrix[0].length; |
| for (var i = 0; i < l1; i++) { |
| for (var j = 0; j < l2; j++) { |
| sum += matrix[i][j] * Math.log(matrix[i][j] + eps); |
| } |
| } |
| return -sum; |
| }; |
| |
| MLStatMatrix.mean = function mean(matrix, dimension) { |
| if (typeof (dimension) === 'undefined') { |
| dimension = 0; |
| } |
| var rows = matrix.length, |
| cols = matrix[0].length, |
| theMean, N, i, j; |
| |
| if (dimension === -1) { |
| theMean = [0]; |
| N = rows * cols; |
| for (i = 0; i < rows; i++) { |
| for (j = 0; j < cols; j++) { |
| theMean[0] += matrix[i][j]; |
| } |
| } |
| theMean[0] /= N; |
| } else if (dimension === 0) { |
| theMean = new Array(cols); |
| N = rows; |
| for (j = 0; j < cols; j++) { |
| theMean[j] = 0; |
| for (i = 0; i < rows; i++) { |
| theMean[j] += matrix[i][j]; |
| } |
| theMean[j] /= N; |
| } |
| } else if (dimension === 1) { |
| theMean = new Array(rows); |
| N = cols; |
| for (j = 0; j < rows; j++) { |
| theMean[j] = 0; |
| for (i = 0; i < cols; i++) { |
| theMean[j] += matrix[j][i]; |
| } |
| theMean[j] /= N; |
| } |
| } else { |
| throw new Error('Invalid dimension'); |
| } |
| return theMean; |
| }; |
| |
| MLStatMatrix.sum = function sum(matrix, dimension) { |
| if (typeof (dimension) === 'undefined') { |
| dimension = 0; |
| } |
| var rows = matrix.length, |
| cols = matrix[0].length, |
| theSum, i, j; |
| |
| if (dimension === -1) { |
| theSum = [0]; |
| for (i = 0; i < rows; i++) { |
| for (j = 0; j < cols; j++) { |
| theSum[0] += matrix[i][j]; |
| } |
| } |
| } else if (dimension === 0) { |
| theSum = new Array(cols); |
| for (j = 0; j < cols; j++) { |
| theSum[j] = 0; |
| for (i = 0; i < rows; i++) { |
| theSum[j] += matrix[i][j]; |
| } |
| } |
| } else if (dimension === 1) { |
| theSum = new Array(rows); |
| for (j = 0; j < rows; j++) { |
| theSum[j] = 0; |
| for (i = 0; i < cols; i++) { |
| theSum[j] += matrix[j][i]; |
| } |
| } |
| } else { |
| throw new Error('Invalid dimension'); |
| } |
| return theSum; |
| }; |
| |
| MLStatMatrix.product = function product(matrix, dimension) { |
| if (typeof (dimension) === 'undefined') { |
| dimension = 0; |
| } |
| var rows = matrix.length, |
| cols = matrix[0].length, |
| theProduct, i, j; |
| |
| if (dimension === -1) { |
| theProduct = [1]; |
| for (i = 0; i < rows; i++) { |
| for (j = 0; j < cols; j++) { |
| theProduct[0] *= matrix[i][j]; |
| } |
| } |
| } else if (dimension === 0) { |
| theProduct = new Array(cols); |
| for (j = 0; j < cols; j++) { |
| theProduct[j] = 1; |
| for (i = 0; i < rows; i++) { |
| theProduct[j] *= matrix[i][j]; |
| } |
| } |
| } else if (dimension === 1) { |
| theProduct = new Array(rows); |
| for (j = 0; j < rows; j++) { |
| theProduct[j] = 1; |
| for (i = 0; i < cols; i++) { |
| theProduct[j] *= matrix[j][i]; |
| } |
| } |
| } else { |
| throw new Error('Invalid dimension'); |
| } |
| return theProduct; |
| }; |
| |
| MLStatMatrix.standardDeviation = function standardDeviation(matrix, means, unbiased) { |
| var vari = MLStatMatrix.variance(matrix, means, unbiased), l = vari.length; |
| for (var i = 0; i < l; i++) { |
| vari[i] = Math.sqrt(vari[i]); |
| } |
| return vari; |
| }; |
| |
| MLStatMatrix.variance = function variance(matrix, means, unbiased) { |
| if (typeof (unbiased) === 'undefined') { |
| unbiased = true; |
| } |
| means = means || MLStatMatrix.mean(matrix); |
| var rows = matrix.length; |
| if (rows === 0) return []; |
| var cols = matrix[0].length; |
| var vari = new Array(cols); |
| |
| for (var j = 0; j < cols; j++) { |
| var sum1 = 0, sum2 = 0, x = 0; |
| for (var i = 0; i < rows; i++) { |
| x = matrix[i][j] - means[j]; |
| sum1 += x; |
| sum2 += x * x; |
| } |
| if (unbiased) { |
| vari[j] = (sum2 - ((sum1 * sum1) / rows)) / (rows - 1); |
| } else { |
| vari[j] = (sum2 - ((sum1 * sum1) / rows)) / rows; |
| } |
| } |
| return vari; |
| }; |
| |
| MLStatMatrix.median = function median(matrix) { |
| var rows = matrix.length, cols = matrix[0].length; |
| var medians = new Array(cols); |
| |
| for (var i = 0; i < cols; i++) { |
| var data = new Array(rows); |
| for (var j = 0; j < rows; j++) { |
| data[j] = matrix[j][i]; |
| } |
| data.sort(compareNumbers); |
| var N = data.length; |
| if (N % 2 === 0) { |
| medians[i] = (data[N / 2] + data[(N / 2) - 1]) * 0.5; |
| } else { |
| medians[i] = data[Math.floor(N / 2)]; |
| } |
| } |
| return medians; |
| }; |
| |
| MLStatMatrix.mode = function mode(matrix) { |
| var rows = matrix.length, |
| cols = matrix[0].length, |
| modes = new Array(cols), |
| i, j; |
| for (i = 0; i < cols; i++) { |
| var itemCount = new Array(rows); |
| for (var k = 0; k < rows; k++) { |
| itemCount[k] = 0; |
| } |
| var itemArray = new Array(rows); |
| var count = 0; |
| |
| for (j = 0; j < rows; j++) { |
| var index = itemArray.indexOf(matrix[j][i]); |
| if (index >= 0) { |
| itemCount[index]++; |
| } else { |
| itemArray[count] = matrix[j][i]; |
| itemCount[count] = 1; |
| count++; |
| } |
| } |
| |
| var maxValue = 0, maxIndex = 0; |
| for (j = 0; j < count; j++) { |
| if (itemCount[j] > maxValue) { |
| maxValue = itemCount[j]; |
| maxIndex = j; |
| } |
| } |
| |
| modes[i] = itemArray[maxIndex]; |
| } |
| return modes; |
| }; |
| |
| MLStatMatrix.skewness = function skewness(matrix, unbiased) { |
| if (typeof (unbiased) === 'undefined') unbiased = true; |
| var means = MLStatMatrix.mean(matrix); |
| var n = matrix.length, l = means.length; |
| var skew = new Array(l); |
| |
| for (var j = 0; j < l; j++) { |
| var s2 = 0, s3 = 0; |
| for (var i = 0; i < n; i++) { |
| var dev = matrix[i][j] - means[j]; |
| s2 += dev * dev; |
| s3 += dev * dev * dev; |
| } |
| |
| var m2 = s2 / n; |
| var m3 = s3 / n; |
| var g = m3 / Math.pow(m2, 3 / 2); |
| |
| if (unbiased) { |
| var a = Math.sqrt(n * (n - 1)); |
| var b = n - 2; |
| skew[j] = (a / b) * g; |
| } else { |
| skew[j] = g; |
| } |
| } |
| return skew; |
| }; |
| |
| MLStatMatrix.kurtosis = function kurtosis(matrix, unbiased) { |
| if (typeof (unbiased) === 'undefined') unbiased = true; |
| var means = MLStatMatrix.mean(matrix); |
| var n = matrix.length, m = matrix[0].length; |
| var kurt = new Array(m); |
| |
| for (var j = 0; j < m; j++) { |
| var s2 = 0, s4 = 0; |
| for (var i = 0; i < n; i++) { |
| var dev = matrix[i][j] - means[j]; |
| s2 += dev * dev; |
| s4 += dev * dev * dev * dev; |
| } |
| var m2 = s2 / n; |
| var m4 = s4 / n; |
| |
| if (unbiased) { |
| var v = s2 / (n - 1); |
| var a = (n * (n + 1)) / ((n - 1) * (n - 2) * (n - 3)); |
| var b = s4 / (v * v); |
| var c = ((n - 1) * (n - 1)) / ((n - 2) * (n - 3)); |
| kurt[j] = a * b - 3 * c; |
| } else { |
| kurt[j] = m4 / (m2 * m2) - 3; |
| } |
| } |
| return kurt; |
| }; |
| |
| MLStatMatrix.standardError = function standardError(matrix) { |
| var samples = matrix.length; |
| var standardDeviations = MLStatMatrix.standardDeviation(matrix); |
| var l = standardDeviations.length; |
| var standardErrors = new Array(l); |
| var sqrtN = Math.sqrt(samples); |
| |
| for (var i = 0; i < l; i++) { |
| standardErrors[i] = standardDeviations[i] / sqrtN; |
| } |
| return standardErrors; |
| }; |
| |
| MLStatMatrix.covariance = function covariance(matrix, dimension) { |
| return MLStatMatrix.scatter(matrix, undefined, dimension); |
| }; |
| |
| MLStatMatrix.scatter = function scatter(matrix, divisor, dimension) { |
| if (typeof (dimension) === 'undefined') { |
| dimension = 0; |
| } |
| if (typeof (divisor) === 'undefined') { |
| if (dimension === 0) { |
| divisor = matrix.length - 1; |
| } else if (dimension === 1) { |
| divisor = matrix[0].length - 1; |
| } |
| } |
| var means = MLStatMatrix.mean(matrix, dimension); |
| var rows = matrix.length; |
| if (rows === 0) { |
| return [[]]; |
| } |
| var cols = matrix[0].length, |
| cov, i, j, s, k; |
| |
| if (dimension === 0) { |
| cov = new Array(cols); |
| for (i = 0; i < cols; i++) { |
| cov[i] = new Array(cols); |
| } |
| for (i = 0; i < cols; i++) { |
| for (j = i; j < cols; j++) { |
| s = 0; |
| for (k = 0; k < rows; k++) { |
| s += (matrix[k][j] - means[j]) * (matrix[k][i] - means[i]); |
| } |
| s /= divisor; |
| cov[i][j] = s; |
| cov[j][i] = s; |
| } |
| } |
| } else if (dimension === 1) { |
| cov = new Array(rows); |
| for (i = 0; i < rows; i++) { |
| cov[i] = new Array(rows); |
| } |
| for (i = 0; i < rows; i++) { |
| for (j = i; j < rows; j++) { |
| s = 0; |
| for (k = 0; k < cols; k++) { |
| s += (matrix[j][k] - means[j]) * (matrix[i][k] - means[i]); |
| } |
| s /= divisor; |
| cov[i][j] = s; |
| cov[j][i] = s; |
| } |
| } |
| } else { |
| throw new Error('Invalid dimension'); |
| } |
| |
| return cov; |
| }; |
| |
| MLStatMatrix.correlation = function correlation(matrix) { |
| var means = MLStatMatrix.mean(matrix), |
| standardDeviations = MLStatMatrix.standardDeviation(matrix, true, means), |
| scores = MLStatMatrix.zScores(matrix, means, standardDeviations), |
| rows = matrix.length, |
| cols = matrix[0].length, |
| i, j; |
| |
| var cor = new Array(cols); |
| for (i = 0; i < cols; i++) { |
| cor[i] = new Array(cols); |
| } |
| for (i = 0; i < cols; i++) { |
| for (j = i; j < cols; j++) { |
| var c = 0; |
| for (var k = 0, l = scores.length; k < l; k++) { |
| c += scores[k][j] * scores[k][i]; |
| } |
| c /= rows - 1; |
| cor[i][j] = c; |
| cor[j][i] = c; |
| } |
| } |
| return cor; |
| }; |
| |
| MLStatMatrix.zScores = function zScores(matrix, means, standardDeviations) { |
| means = means || MLStatMatrix.mean(matrix); |
| if (typeof (standardDeviations) === 'undefined') standardDeviations = MLStatMatrix.standardDeviation(matrix, true, means); |
| return MLStatMatrix.standardize(MLStatMatrix.center(matrix, means, false), standardDeviations, true); |
| }; |
| |
| MLStatMatrix.center = function center(matrix, means, inPlace) { |
| means = means || MLStatMatrix.mean(matrix); |
| var result = matrix, |
| l = matrix.length, |
| i, j, jj; |
| |
| if (!inPlace) { |
| result = new Array(l); |
| for (i = 0; i < l; i++) { |
| result[i] = new Array(matrix[i].length); |
| } |
| } |
| |
| for (i = 0; i < l; i++) { |
| var row = result[i]; |
| for (j = 0, jj = row.length; j < jj; j++) { |
| row[j] = matrix[i][j] - means[j]; |
| } |
| } |
| return result; |
| }; |
| |
| MLStatMatrix.standardize = function standardize(matrix, standardDeviations, inPlace) { |
| if (typeof (standardDeviations) === 'undefined') standardDeviations = MLStatMatrix.standardDeviation(matrix); |
| var result = matrix, |
| l = matrix.length, |
| i, j, jj; |
| |
| if (!inPlace) { |
| result = new Array(l); |
| for (i = 0; i < l; i++) { |
| result[i] = new Array(matrix[i].length); |
| } |
| } |
| |
| for (i = 0; i < l; i++) { |
| var resultRow = result[i]; |
| var sourceRow = matrix[i]; |
| for (j = 0, jj = resultRow.length; j < jj; j++) { |
| if (standardDeviations[j] !== 0 && !isNaN(standardDeviations[j])) { |
| resultRow[j] = sourceRow[j] / standardDeviations[j]; |
| } |
| } |
| } |
| return result; |
| }; |
| |
| MLStatMatrix.weightedVariance = function weightedVariance(matrix, weights) { |
| var means = MLStatMatrix.mean(matrix); |
| var rows = matrix.length; |
| if (rows === 0) return []; |
| var cols = matrix[0].length; |
| var vari = new Array(cols); |
| |
| for (var j = 0; j < cols; j++) { |
| var sum = 0; |
| var a = 0, b = 0; |
| |
| for (var i = 0; i < rows; i++) { |
| var z = matrix[i][j] - means[j]; |
| var w = weights[i]; |
| |
| sum += w * (z * z); |
| b += w; |
| a += w * w; |
| } |
| |
| vari[j] = sum * (b / (b * b - a)); |
| } |
| |
| return vari; |
| }; |
| |
| MLStatMatrix.weightedMean = function weightedMean(matrix, weights, dimension) { |
| if (typeof (dimension) === 'undefined') { |
| dimension = 0; |
| } |
| var rows = matrix.length; |
| if (rows === 0) return []; |
| var cols = matrix[0].length, |
| means, i, ii, j, w, row; |
| |
| if (dimension === 0) { |
| means = new Array(cols); |
| for (i = 0; i < cols; i++) { |
| means[i] = 0; |
| } |
| for (i = 0; i < rows; i++) { |
| row = matrix[i]; |
| w = weights[i]; |
| for (j = 0; j < cols; j++) { |
| means[j] += row[j] * w; |
| } |
| } |
| } else if (dimension === 1) { |
| means = new Array(rows); |
| for (i = 0; i < rows; i++) { |
| means[i] = 0; |
| } |
| for (j = 0; j < rows; j++) { |
| row = matrix[j]; |
| w = weights[j]; |
| for (i = 0; i < cols; i++) { |
| means[j] += row[i] * w; |
| } |
| } |
| } else { |
| throw new Error('Invalid dimension'); |
| } |
| |
| var weightSum = arrayStat.sum(weights); |
| if (weightSum !== 0) { |
| for (i = 0, ii = means.length; i < ii; i++) { |
| means[i] /= weightSum; |
| } |
| } |
| return means; |
| }; |
| |
| MLStatMatrix.weightedCovariance = function weightedCovariance(matrix, weights, means, dimension) { |
| dimension = dimension || 0; |
| means = means || MLStatMatrix.weightedMean(matrix, weights, dimension); |
| var s1 = 0, s2 = 0; |
| for (var i = 0, ii = weights.length; i < ii; i++) { |
| s1 += weights[i]; |
| s2 += weights[i] * weights[i]; |
| } |
| var factor = s1 / (s1 * s1 - s2); |
| return MLStatMatrix.weightedScatter(matrix, weights, means, factor, dimension); |
| }; |
| |
| MLStatMatrix.weightedScatter = function weightedScatter(matrix, weights, means, factor, dimension) { |
| dimension = dimension || 0; |
| means = means || MLStatMatrix.weightedMean(matrix, weights, dimension); |
| if (typeof (factor) === 'undefined') { |
| factor = 1; |
| } |
| var rows = matrix.length; |
| if (rows === 0) { |
| return [[]]; |
| } |
| var cols = matrix[0].length, |
| cov, i, j, k, s; |
| |
| if (dimension === 0) { |
| cov = new Array(cols); |
| for (i = 0; i < cols; i++) { |
| cov[i] = new Array(cols); |
| } |
| for (i = 0; i < cols; i++) { |
| for (j = i; j < cols; j++) { |
| s = 0; |
| for (k = 0; k < rows; k++) { |
| s += weights[k] * (matrix[k][j] - means[j]) * (matrix[k][i] - means[i]); |
| } |
| cov[i][j] = s * factor; |
| cov[j][i] = s * factor; |
| } |
| } |
| } else if (dimension === 1) { |
| cov = new Array(rows); |
| for (i = 0; i < rows; i++) { |
| cov[i] = new Array(rows); |
| } |
| for (i = 0; i < rows; i++) { |
| for (j = i; j < rows; j++) { |
| s = 0; |
| for (k = 0; k < cols; k++) { |
| s += weights[k] * (matrix[j][k] - means[j]) * (matrix[i][k] - means[i]); |
| } |
| cov[i][j] = s * factor; |
| cov[j][i] = s * factor; |
| } |
| } |
| } else { |
| throw new Error('Invalid dimension'); |
| } |
| |
| return cov; |
| }; |
| } |
| |
| // ml-stat index.js |
| const MLStat = {}; |
| { |
| MLStat.array = MLStatArray; |
| MLStat.matrix = MLStatMatrix; |
| } |
| |
| |
| // ml-array-utils ArrayUtils.js |
| const MLArrayUtilsArrayUtils = {}; |
| { |
| const Stat = MLStat.array; |
| /** |
| * Function that returns an array of points given 1D array as follows: |
| * |
| * [x1, y1, .. , x2, y2, ..] |
| * |
| * And receive the number of dimensions of each point. |
| * @param array |
| * @param dimensions |
| * @returns {Array} - Array of points. |
| */ |
| function coordArrayToPoints(array, dimensions) { |
| if(array.length % dimensions !== 0) { |
| throw new RangeError('Dimensions number must be accordance with the size of the array.'); |
| } |
| |
| var length = array.length / dimensions; |
| var pointsArr = new Array(length); |
| |
| var k = 0; |
| for(var i = 0; i < array.length; i += dimensions) { |
| var point = new Array(dimensions); |
| for(var j = 0; j < dimensions; ++j) { |
| point[j] = array[i + j]; |
| } |
| |
| pointsArr[k] = point; |
| k++; |
| } |
| |
| return pointsArr; |
| } |
| |
| |
| /** |
| * Function that given an array as follows: |
| * [x1, y1, .. , x2, y2, ..] |
| * |
| * Returns an array as follows: |
| * [[x1, x2, ..], [y1, y2, ..], [ .. ]] |
| * |
| * And receives the number of dimensions of each coordinate. |
| * @param array |
| * @param dimensions |
| * @returns {Array} - Matrix of coordinates |
| */ |
| function coordArrayToCoordMatrix(array, dimensions) { |
| if(array.length % dimensions !== 0) { |
| throw new RangeError('Dimensions number must be accordance with the size of the array.'); |
| } |
| |
| var coordinatesArray = new Array(dimensions); |
| var points = array.length / dimensions; |
| for (var i = 0; i < coordinatesArray.length; i++) { |
| coordinatesArray[i] = new Array(points); |
| } |
| |
| for(i = 0; i < array.length; i += dimensions) { |
| for(var j = 0; j < dimensions; ++j) { |
| var currentPoint = Math.floor(i / dimensions); |
| coordinatesArray[j][currentPoint] = array[i + j]; |
| } |
| } |
| |
| return coordinatesArray; |
| } |
| |
| /** |
| * Function that receives a coordinate matrix as follows: |
| * [[x1, x2, ..], [y1, y2, ..], [ .. ]] |
| * |
| * Returns an array of coordinates as follows: |
| * [x1, y1, .. , x2, y2, ..] |
| * |
| * @param coordMatrix |
| * @returns {Array} |
| */ |
| function coordMatrixToCoordArray(coordMatrix) { |
| var coodinatesArray = new Array(coordMatrix.length * coordMatrix[0].length); |
| var k = 0; |
| for(var i = 0; i < coordMatrix[0].length; ++i) { |
| for(var j = 0; j < coordMatrix.length; ++j) { |
| coodinatesArray[k] = coordMatrix[j][i]; |
| ++k; |
| } |
| } |
| |
| return coodinatesArray; |
| } |
| |
| /** |
| * Tranpose a matrix, this method is for coordMatrixToPoints and |
| * pointsToCoordMatrix, that because only transposing the matrix |
| * you can change your representation. |
| * |
| * @param matrix |
| * @returns {Array} |
| */ |
| function transpose(matrix) { |
| var resultMatrix = new Array(matrix[0].length); |
| for(var i = 0; i < resultMatrix.length; ++i) { |
| resultMatrix[i] = new Array(matrix.length); |
| } |
| |
| for (i = 0; i < matrix.length; ++i) { |
| for(var j = 0; j < matrix[0].length; ++j) { |
| resultMatrix[j][i] = matrix[i][j]; |
| } |
| } |
| |
| return resultMatrix; |
| } |
| |
| /** |
| * Function that transform an array of points into a coordinates array |
| * as follows: |
| * [x1, y1, .. , x2, y2, ..] |
| * |
| * @param points |
| * @returns {Array} |
| */ |
| function pointsToCoordArray(points) { |
| var coodinatesArray = new Array(points.length * points[0].length); |
| var k = 0; |
| for(var i = 0; i < points.length; ++i) { |
| for(var j = 0; j < points[0].length; ++j) { |
| coodinatesArray[k] = points[i][j]; |
| ++k; |
| } |
| } |
| |
| return coodinatesArray; |
| } |
| |
| /** |
| * Apply the dot product between the smaller vector and a subsets of the |
| * largest one. |
| * |
| * @param firstVector |
| * @param secondVector |
| * @returns {Array} each dot product of size of the difference between the |
| * larger and the smallest one. |
| */ |
| function applyDotProduct(firstVector, secondVector) { |
| var largestVector, smallestVector; |
| if(firstVector.length <= secondVector.length) { |
| smallestVector = firstVector; |
| largestVector = secondVector; |
| } else { |
| smallestVector = secondVector; |
| largestVector = firstVector; |
| } |
| |
| var difference = largestVector.length - smallestVector.length + 1; |
| var dotProductApplied = new Array(difference); |
| |
| for (var i = 0; i < difference; ++i) { |
| var sum = 0; |
| for (var j = 0; j < smallestVector.length; ++j) { |
| sum += smallestVector[j] * largestVector[i + j]; |
| } |
| dotProductApplied[i] = sum; |
| } |
| |
| return dotProductApplied; |
| } |
| /** |
| * To scale the input array between the specified min and max values. The operation is performed inplace |
| * if the options.inplace is specified. If only one of the min or max parameters is specified, then the scaling |
| * will multiply the input array by min/min(input) or max/max(input) |
| * @param input |
| * @param options |
| * @returns {*} |
| */ |
| function scale(input, options){ |
| var y; |
| if(options.inPlace){ |
| y = input; |
| } |
| else{ |
| y = new Array(input.length); |
| } |
| const max = options.max; |
| const min = options.min; |
| if(typeof max === "number"){ |
| if(typeof min === "number"){ |
| var minMax = Stat.minMax(input); |
| var factor = (max - min)/(minMax.max-minMax.min); |
| for(var i=0;i< y.length;i++){ |
| y[i]=(input[i]-minMax.min)*factor+min; |
| } |
| } |
| else{ |
| var currentMin = Stat.max(input); |
| var factor = max/currentMin; |
| for(var i=0;i< y.length;i++){ |
| y[i] = input[i]*factor; |
| } |
| } |
| } |
| else{ |
| if(typeof min === "number"){ |
| var currentMin = Stat.min(input); |
| var factor = min/currentMin; |
| for(var i=0;i< y.length;i++){ |
| y[i] = input[i]*factor; |
| } |
| } |
| } |
| return y; |
| } |
| |
| MLArrayUtilsArrayUtils.coordArrayToPoints = coordArrayToPoints; |
| MLArrayUtilsArrayUtils.coordArrayToCoordMatrix = coordArrayToCoordMatrix; |
| MLArrayUtilsArrayUtils.coordMatrixToCoordArray = coordMatrixToCoordArray; |
| MLArrayUtilsArrayUtils.coordMatrixToPoints = transpose; |
| MLArrayUtilsArrayUtils.pointsToCoordArray = pointsToCoordArray; |
| MLArrayUtilsArrayUtils.pointsToCoordMatrix = transpose; |
| MLArrayUtilsArrayUtils.applyDotProduct = applyDotProduct; |
| MLArrayUtilsArrayUtils.scale = scale; |
| } |
| |
| |
| // ml-array-utils getEquallySpaced.js |
| const MLArrayUtilsGetEquallySpaced = {}; |
| { |
| /** |
| * |
| * Function that returns a Number array of equally spaced numberOfPoints |
| * containing a representation of intensities of the spectra arguments x |
| * and y. |
| * |
| * The options parameter contains an object in the following form: |
| * from: starting point |
| * to: last point |
| * numberOfPoints: number of points between from and to |
| * variant: "slot" or "smooth" - smooth is the default option |
| * |
| * The slot variant consist that each point in the new array is calculated |
| * averaging the existing points between the slot that belongs to the current |
| * value. The smooth variant is the same but takes the integral of the range |
| * of the slot and divide by the step size between two points in the new array. |
| * |
| * @param x - sorted increasing x values |
| * @param y |
| * @param options |
| * @returns {Array} new array with the equally spaced data. |
| * |
| */ |
| function getEquallySpacedData(x, y, options) { |
| if (x.length>1 && x[0]>x[1]) { |
| x=x.slice().reverse(); |
| y=y.slice().reverse(); |
| } |
| |
| var xLength = x.length; |
| if(xLength !== y.length) |
| throw new RangeError("the x and y vector doesn't have the same size."); |
| |
| if (options === undefined) options = {}; |
| |
| var from = options.from === undefined ? x[0] : options.from |
| if (isNaN(from) || !isFinite(from)) { |
| throw new RangeError("'From' value must be a number"); |
| } |
| var to = options.to === undefined ? x[x.length - 1] : options.to; |
| if (isNaN(to) || !isFinite(to)) { |
| throw new RangeError("'To' value must be a number"); |
| } |
| |
| var reverse = from > to; |
| if(reverse) { |
| var temp = from; |
| from = to; |
| to = temp; |
| } |
| |
| var numberOfPoints = options.numberOfPoints === undefined ? 100 : options.numberOfPoints; |
| if (isNaN(numberOfPoints) || !isFinite(numberOfPoints)) { |
| throw new RangeError("'Number of points' value must be a number"); |
| } |
| if(numberOfPoints < 1) |
| throw new RangeError("the number of point must be higher than 1"); |
| |
| var algorithm = options.variant === "slot" ? "slot" : "smooth"; // default value: smooth |
| |
| var output = algorithm === "slot" ? getEquallySpacedSlot(x, y, from, to, numberOfPoints) : getEquallySpacedSmooth(x, y, from, to, numberOfPoints); |
| |
| return reverse ? output.reverse() : output; |
| } |
| |
| /** |
| * function that retrieves the getEquallySpacedData with the variant "smooth" |
| * |
| * @param x |
| * @param y |
| * @param from - Initial point |
| * @param to - Final point |
| * @param numberOfPoints |
| * @returns {Array} - Array of y's equally spaced with the variant "smooth" |
| */ |
| function getEquallySpacedSmooth(x, y, from, to, numberOfPoints) { |
| var xLength = x.length; |
| |
| var step = (to - from) / (numberOfPoints - 1); |
| var halfStep = step / 2; |
| |
| var start = from - halfStep; |
| var output = new Array(numberOfPoints); |
| |
| var initialOriginalStep = x[1] - x[0]; |
| var lastOriginalStep = x[x.length - 1] - x[x.length - 2]; |
| |
| // Init main variables |
| var min = start; |
| var max = start + step; |
| |
| var previousX = Number.MIN_VALUE; |
| var previousY = 0; |
| var nextX = x[0] - initialOriginalStep; |
| var nextY = 0; |
| |
| var currentValue = 0; |
| var slope = 0; |
| var intercept = 0; |
| var sumAtMin = 0; |
| var sumAtMax = 0; |
| |
| var i = 0; // index of input |
| var j = 0; // index of output |
| |
| function getSlope(x0, y0, x1, y1) { |
| return (y1 - y0) / (x1 - x0); |
| } |
| |
| main: while(true) { |
| while (nextX - max >= 0) { |
| // no overlap with original point, just consume current value |
| var add = integral(0, max - previousX, slope, previousY); |
| sumAtMax = currentValue + add; |
| |
| output[j] = (sumAtMax - sumAtMin) / step; |
| j++; |
| |
| if (j === numberOfPoints) |
| break main; |
| |
| min = max; |
| max += step; |
| sumAtMin = sumAtMax; |
| } |
| |
| if(previousX <= min && min <= nextX) { |
| add = integral(0, min - previousX, slope, previousY); |
| sumAtMin = currentValue + add; |
| } |
| |
| currentValue += integral(previousX, nextX, slope, intercept); |
| |
| previousX = nextX; |
| previousY = nextY; |
| |
| if (i < xLength) { |
| nextX = x[i]; |
| nextY = y[i]; |
| i++; |
| } else if (i === xLength) { |
| nextX += lastOriginalStep; |
| nextY = 0; |
| } |
| // updating parameters |
| slope = getSlope(previousX, previousY, nextX, nextY); |
| intercept = -slope*previousX + previousY; |
| } |
| |
| return output; |
| } |
| |
| /** |
| * function that retrieves the getEquallySpacedData with the variant "slot" |
| * |
| * @param x |
| * @param y |
| * @param from - Initial point |
| * @param to - Final point |
| * @param numberOfPoints |
| * @returns {Array} - Array of y's equally spaced with the variant "slot" |
| */ |
| function getEquallySpacedSlot(x, y, from, to, numberOfPoints) { |
| var xLength = x.length; |
| |
| var step = (to - from) / (numberOfPoints - 1); |
| var halfStep = step / 2; |
| var lastStep = x[x.length - 1] - x[x.length - 2]; |
| |
| var start = from - halfStep; |
| var output = new Array(numberOfPoints); |
| |
| // Init main variables |
| var min = start; |
| var max = start + step; |
| |
| var previousX = -Number.MAX_VALUE; |
| var previousY = 0; |
| var nextX = x[0]; |
| var nextY = y[0]; |
| var frontOutsideSpectra = 0; |
| var backOutsideSpectra = true; |
| |
| var currentValue = 0; |
| |
| // for slot algorithm |
| var currentPoints = 0; |
| |
| var i = 1; // index of input |
| var j = 0; // index of output |
| |
| main: while(true) { |
| if (previousX>=nextX) throw (new Error('x must be an increasing serie')); |
| while (previousX - max > 0) { |
| // no overlap with original point, just consume current value |
| if(backOutsideSpectra) { |
| currentPoints++; |
| backOutsideSpectra = false; |
| } |
| |
| output[j] = currentPoints <= 0 ? 0 : currentValue / currentPoints; |
| j++; |
| |
| if (j === numberOfPoints) |
| break main; |
| |
| min = max; |
| max += step; |
| currentValue = 0; |
| currentPoints = 0; |
| } |
| |
| if(previousX > min) { |
| currentValue += previousY; |
| currentPoints++; |
| } |
| |
| if(previousX === -Number.MAX_VALUE || frontOutsideSpectra > 1) |
| currentPoints--; |
| |
| previousX = nextX; |
| previousY = nextY; |
| |
| if (i < xLength) { |
| nextX = x[i]; |
| nextY = y[i]; |
| i++; |
| } else { |
| nextX += lastStep; |
| nextY = 0; |
| frontOutsideSpectra++; |
| } |
| } |
| |
| return output; |
| } |
| /** |
| * Function that calculates the integral of the line between two |
| * x-coordinates, given the slope and intercept of the line. |
| * |
| * @param x0 |
| * @param x1 |
| * @param slope |
| * @param intercept |
| * @returns {number} integral value. |
| */ |
| function integral(x0, x1, slope, intercept) { |
| return (0.5 * slope * x1 * x1 + intercept * x1) - (0.5 * slope * x0 * x0 + intercept * x0); |
| } |
| |
| MLArrayUtilsGetEquallySpaced.getEquallySpacedData = getEquallySpacedData; |
| MLArrayUtilsGetEquallySpaced.integral = integral; |
| } |
| |
| |
| // ml-array-utils snv.js |
| const MLArrayUtilsSNV = {}; |
| { |
| MLArrayUtilsSNV.SNV = SNV; |
| let Stat = MLStat.array; |
| |
| /** |
| * Function that applies the standard normal variate (SNV) to an array of values. |
| * |
| * @param data - Array of values. |
| * @returns {Array} - applied the SNV. |
| */ |
| function SNV(data) { |
| var mean = Stat.mean(data); |
| var std = Stat.standardDeviation(data); |
| var result = data.slice(); |
| for (var i = 0; i < data.length; i++) { |
| result[i] = (result[i] - mean) / std; |
| } |
| return result; |
| } |
| } |
| |
| // ml-array-utils index.js |
| const MLArrayUtils = {}; |
| { |
| MLArrayUtils.getEquallySpacedData = MLArrayUtilsGetEquallySpaced.getEquallySpacedData; |
| MLArrayUtils.SNV = MLArrayUtilsSNV.SNV; |
| } |
| |
| |
| |
| // do this early so things can use it. This is from ml-matrix src/matrix.js |
| const MLMatrixMatrix = {}; |
| |
| // ml-matrix src/util.js |
| const MLMatrixUtil = {}; |
| { |
| let exports = MLMatrixUtil; |
| let Matrix = MLMatrixMatrix; |
| |
| /** |
| * @private |
| * Check that a row index is not out of bounds |
| * @param {Matrix} matrix |
| * @param {number} index |
| * @param {boolean} [outer] |
| */ |
| exports.checkRowIndex = function checkRowIndex(matrix, index, outer) { |
| var max = outer ? matrix.rows : matrix.rows - 1; |
| if (index < 0 || index > max) { |
| throw new RangeError('Row index out of range'); |
| } |
| }; |
| |
| /** |
| * @private |
| * Check that a column index is not out of bounds |
| * @param {Matrix} matrix |
| * @param {number} index |
| * @param {boolean} [outer] |
| */ |
| exports.checkColumnIndex = function checkColumnIndex(matrix, index, outer) { |
| var max = outer ? matrix.columns : matrix.columns - 1; |
| if (index < 0 || index > max) { |
| throw new RangeError('Column index out of range'); |
| } |
| }; |
| |
| /** |
| * @private |
| * Check that the provided vector is an array with the right length |
| * @param {Matrix} matrix |
| * @param {Array|Matrix} vector |
| * @return {Array} |
| * @throws {RangeError} |
| */ |
| exports.checkRowVector = function checkRowVector(matrix, vector) { |
| if (vector.to1DArray) { |
| vector = vector.to1DArray(); |
| } |
| if (vector.length !== matrix.columns) { |
| throw new RangeError('vector size must be the same as the number of columns'); |
| } |
| return vector; |
| }; |
| |
| /** |
| * @private |
| * Check that the provided vector is an array with the right length |
| * @param {Matrix} matrix |
| * @param {Array|Matrix} vector |
| * @return {Array} |
| * @throws {RangeError} |
| */ |
| exports.checkColumnVector = function checkColumnVector(matrix, vector) { |
| if (vector.to1DArray) { |
| vector = vector.to1DArray(); |
| } |
| if (vector.length !== matrix.rows) { |
| throw new RangeError('vector size must be the same as the number of rows'); |
| } |
| return vector; |
| }; |
| |
| exports.checkIndices = function checkIndices(matrix, rowIndices, columnIndices) { |
| var rowOut = rowIndices.some(r => { |
| return r < 0 || r >= matrix.rows; |
| |
| }); |
| |
| var columnOut = columnIndices.some(c => { |
| return c < 0 || c >= matrix.columns; |
| }); |
| |
| if (rowOut || columnOut) { |
| throw new RangeError('Indices are out of range'); |
| } |
| |
| if (typeof rowIndices !== 'object' || typeof columnIndices !== 'object') { |
| throw new TypeError('Unexpected type for row/column indices'); |
| } |
| if (!Array.isArray(rowIndices)) rowIndices = Array.from(rowIndices); |
| if (!Array.isArray(columnIndices)) rowIndices = Array.from(columnIndices); |
| |
| return { |
| row: rowIndices, |
| column: columnIndices |
| }; |
| }; |
| |
| exports.checkRange = function checkRange(matrix, startRow, endRow, startColumn, endColumn) { |
| if (arguments.length !== 5) throw new TypeError('Invalid argument type'); |
| var notAllNumbers = Array.from(arguments).slice(1).some(function (arg) { |
| return typeof arg !== 'number'; |
| }); |
| if (notAllNumbers) throw new TypeError('Invalid argument type'); |
| if (startRow > endRow || startColumn > endColumn || startRow < 0 || startRow >= matrix.rows || endRow < 0 || endRow >= matrix.rows || startColumn < 0 || startColumn >= matrix.columns || endColumn < 0 || endColumn >= matrix.columns) { |
| throw new RangeError('Submatrix indices are out of range'); |
| } |
| }; |
| |
| exports.getRange = function getRange(from, to) { |
| var arr = new Array(to - from + 1); |
| for (var i = 0; i < arr.length; i++) { |
| arr[i] = from + i; |
| } |
| return arr; |
| }; |
| |
| exports.sumByRow = function sumByRow(matrix) { |
| var sum = Matrix.Matrix.zeros(matrix.rows, 1); |
| for (var i = 0; i < matrix.rows; ++i) { |
| for (var j = 0; j < matrix.columns; ++j) { |
| sum.set(i, 0, sum.get(i, 0) + matrix.get(i, j)); |
| } |
| } |
| return sum; |
| }; |
| |
| exports.sumByColumn = function sumByColumn(matrix) { |
| var sum = Matrix.Matrix.zeros(1, matrix.columns); |
| for (var i = 0; i < matrix.rows; ++i) { |
| for (var j = 0; j < matrix.columns; ++j) { |
| sum.set(0, j, sum.get(0, j) + matrix.get(i, j)); |
| } |
| } |
| return sum; |
| }; |
| |
| exports.sumAll = function sumAll(matrix) { |
| var v = 0; |
| for (var i = 0; i < matrix.rows; i++) { |
| for (var j = 0; j < matrix.columns; j++) { |
| v += matrix.get(i, j); |
| } |
| } |
| return v; |
| }; |
| } |
| |
| // ml-matrix symbolsspecies.js |
| if (!Symbol.species) { |
| Symbol.species = Symbol.for('@@species'); |
| } |
| |
| |
| // ml-matrix src/dc/util.js |
| const MLMatrixDCUtil = {}; |
| { |
| let exports = MLMatrixDCUtil; |
| exports.hypotenuse = function hypotenuse(a, b) { |
| var r; |
| if (Math.abs(a) > Math.abs(b)) { |
| r = b / a; |
| return Math.abs(a) * Math.sqrt(1 + r * r); |
| } |
| if (b !== 0) { |
| r = a / b; |
| return Math.abs(b) * Math.sqrt(1 + r * r); |
| } |
| return 0; |
| }; |
| |
| // For use in the decomposition algorithms. With big matrices, access time is |
| // too long on elements from array subclass |
| // todo check when it is fixed in v8 |
| // http://jsperf.com/access-and-write-array-subclass |
| exports.getEmpty2DArray = function (rows, columns) { |
| var array = new Array(rows); |
| for (var i = 0; i < rows; i++) { |
| array[i] = new Array(columns); |
| } |
| return array; |
| }; |
| |
| exports.getFilled2DArray = function (rows, columns, value) { |
| var array = new Array(rows); |
| for (var i = 0; i < rows; i++) { |
| array[i] = new Array(columns); |
| for (var j = 0; j < columns; j++) { |
| array[i][j] = value; |
| } |
| } |
| return array; |
| }; |
| } |
| |
| // ml-matrix src/dc/lu.js |
| let MLMatrixDCLU = {}; |
| { |
| let Matrix = MLMatrixMatrix; |
| |
| // https://github.com/lutzroeder/Mapack/blob/master/Source/LuDecomposition.cs |
| function LuDecomposition(matrix) { |
| if (!(this instanceof LuDecomposition)) { |
| return new LuDecomposition(matrix); |
| } |
| |
| matrix = Matrix.Matrix.checkMatrix(matrix); |
| |
| var lu = matrix.clone(), |
| rows = lu.rows, |
| columns = lu.columns, |
| pivotVector = new Array(rows), |
| pivotSign = 1, |
| i, j, k, p, s, t, v, |
| LUrowi, LUcolj, kmax; |
| |
| for (i = 0; i < rows; i++) { |
| pivotVector[i] = i; |
| } |
| |
| LUcolj = new Array(rows); |
| |
| for (j = 0; j < columns; j++) { |
| |
| for (i = 0; i < rows; i++) { |
| LUcolj[i] = lu[i][j]; |
| } |
| |
| for (i = 0; i < rows; i++) { |
| LUrowi = lu[i]; |
| kmax = Math.min(i, j); |
| s = 0; |
| for (k = 0; k < kmax; k++) { |
| s += LUrowi[k] * LUcolj[k]; |
| } |
| LUrowi[j] = LUcolj[i] -= s; |
| } |
| |
| p = j; |
| for (i = j + 1; i < rows; i++) { |
| if (Math.abs(LUcolj[i]) > Math.abs(LUcolj[p])) { |
| p = i; |
| } |
| } |
| |
| if (p !== j) { |
| for (k = 0; k < columns; k++) { |
| t = lu[p][k]; |
| lu[p][k] = lu[j][k]; |
| lu[j][k] = t; |
| } |
| |
| v = pivotVector[p]; |
| pivotVector[p] = pivotVector[j]; |
| pivotVector[j] = v; |
| |
| pivotSign = -pivotSign; |
| } |
| |
| if (j < rows && lu[j][j] !== 0) { |
| for (i = j + 1; i < rows; i++) { |
| lu[i][j] /= lu[j][j]; |
| } |
| } |
| } |
| |
| this.LU = lu; |
| this.pivotVector = pivotVector; |
| this.pivotSign = pivotSign; |
| } |
| |
| LuDecomposition.prototype = { |
| isSingular: function () { |
| var data = this.LU, |
| col = data.columns; |
| for (var j = 0; j < col; j++) { |
| if (data[j][j] === 0) { |
| return true; |
| } |
| } |
| return false; |
| }, |
| get determinant() { |
| var data = this.LU; |
| if (!data.isSquare()) { |
| throw new Error('Matrix must be square'); |
| } |
| var determinant = this.pivotSign, col = data.columns; |
| for (var j = 0; j < col; j++) { |
| determinant *= data[j][j]; |
| } |
| return determinant; |
| }, |
| get lowerTriangularMatrix() { |
| var data = this.LU, |
| rows = data.rows, |
| columns = data.columns, |
| X = new Matrix.Matrix(rows, columns); |
| for (var i = 0; i < rows; i++) { |
| for (var j = 0; j < columns; j++) { |
| if (i > j) { |
| X[i][j] = data[i][j]; |
| } else if (i === j) { |
| X[i][j] = 1; |
| } else { |
| X[i][j] = 0; |
| } |
| } |
| } |
| return X; |
| }, |
| get upperTriangularMatrix() { |
| var data = this.LU, |
| rows = data.rows, |
| columns = data.columns, |
| X = new Matrix.Matrix(rows, columns); |
| for (var i = 0; i < rows; i++) { |
| for (var j = 0; j < columns; j++) { |
| if (i <= j) { |
| X[i][j] = data[i][j]; |
| } else { |
| X[i][j] = 0; |
| } |
| } |
| } |
| return X; |
| }, |
| get pivotPermutationVector() { |
| return this.pivotVector.slice(); |
| }, |
| solve: function (value) { |
| value = Matrix.Matrix.checkMatrix(value); |
| |
| var lu = this.LU, |
| rows = lu.rows; |
| |
| if (rows !== value.rows) { |
| throw new Error('Invalid matrix dimensions'); |
| } |
| if (this.isSingular()) { |
| throw new Error('LU matrix is singular'); |
| } |
| |
| var count = value.columns; |
| var X = value.subMatrixRow(this.pivotVector, 0, count - 1); |
| var columns = lu.columns; |
| var i, j, k; |
| |
| for (k = 0; k < columns; k++) { |
| for (i = k + 1; i < columns; i++) { |
| for (j = 0; j < count; j++) { |
| X[i][j] -= X[k][j] * lu[i][k]; |
| } |
| } |
| } |
| for (k = columns - 1; k >= 0; k--) { |
| for (j = 0; j < count; j++) { |
| X[k][j] /= lu[k][k]; |
| } |
| for (i = 0; i < k; i++) { |
| for (j = 0; j < count; j++) { |
| X[i][j] -= X[k][j] * lu[i][k]; |
| } |
| } |
| } |
| return X; |
| } |
| }; |
| |
| MLMatrixDCLU = LuDecomposition; |
| } |
| |
| |
| // ml-matrix src/dc/svd.js |
| let MLMatrixDCSVD = {}; |
| { |
| let Matrix = MLMatrixMatrix; |
| let util = MLMatrixDCUtil; |
| let hypotenuse = util.hypotenuse; |
| let getFilled2DArray = util.getFilled2DArray; |
| |
| // https://github.com/lutzroeder/Mapack/blob/master/Source/SingularValueDecomposition.cs |
| function SingularValueDecomposition(value, options) { |
| if (!(this instanceof SingularValueDecomposition)) { |
| return new SingularValueDecomposition(value, options); |
| } |
| value = Matrix.Matrix.checkMatrix(value); |
| |
| options = options || {}; |
| |
| var m = value.rows, |
| n = value.columns, |
| nu = Math.min(m, n); |
| |
| var wantu = true, wantv = true; |
| if (options.computeLeftSingularVectors === false) wantu = false; |
| if (options.computeRightSingularVectors === false) wantv = false; |
| var autoTranspose = options.autoTranspose === true; |
| |
| var swapped = false; |
| var a; |
| if (m < n) { |
| if (!autoTranspose) { |
| a = value.clone(); |
| // eslint-disable-next-line no-console |
| console.warn('Computing SVD on a matrix with more columns than rows. Consider enabling autoTranspose'); |
| } else { |
| a = value.transpose(); |
| m = a.rows; |
| n = a.columns; |
| swapped = true; |
| var aux = wantu; |
| wantu = wantv; |
| wantv = aux; |
| } |
| } else { |
| a = value.clone(); |
| } |
| |
| var s = new Array(Math.min(m + 1, n)), |
| U = getFilled2DArray(m, nu, 0), |
| V = getFilled2DArray(n, n, 0), |
| e = new Array(n), |
| work = new Array(m); |
| |
| var nct = Math.min(m - 1, n); |
| var nrt = Math.max(0, Math.min(n - 2, m)); |
| |
| var i, j, k, p, t, ks, f, cs, sn, max, kase, |
| scale, sp, spm1, epm1, sk, ek, b, c, shift, g; |
| |
| for (k = 0, max = Math.max(nct, nrt); k < max; k++) { |
| if (k < nct) { |
| s[k] = 0; |
| for (i = k; i < m; i++) { |
| s[k] = hypotenuse(s[k], a[i][k]); |
| } |
| if (s[k] !== 0) { |
| if (a[k][k] < 0) { |
| s[k] = -s[k]; |
| } |
| for (i = k; i < m; i++) { |
| a[i][k] /= s[k]; |
| } |
| a[k][k] += 1; |
| } |
| s[k] = -s[k]; |
| } |
| |
| for (j = k + 1; j < n; j++) { |
| if ((k < nct) && (s[k] !== 0)) { |
| t = 0; |
| for (i = k; i < m; i++) { |
| t += a[i][k] * a[i][j]; |
| } |
| t = -t / a[k][k]; |
| for (i = k; i < m; i++) { |
| a[i][j] += t * a[i][k]; |
| } |
| } |
| e[j] = a[k][j]; |
| } |
| |
| if (wantu && (k < nct)) { |
| for (i = k; i < m; i++) { |
| U[i][k] = a[i][k]; |
| } |
| } |
| |
| if (k < nrt) { |
| e[k] = 0; |
| for (i = k + 1; i < n; i++) { |
| e[k] = hypotenuse(e[k], e[i]); |
| } |
| if (e[k] !== 0) { |
| if (e[k + 1] < 0) { |
| e[k] = 0 - e[k]; |
| } |
| for (i = k + 1; i < n; i++) { |
| e[i] /= e[k]; |
| } |
| e[k + 1] += 1; |
| } |
| e[k] = -e[k]; |
| if ((k + 1 < m) && (e[k] !== 0)) { |
| for (i = k + 1; i < m; i++) { |
| work[i] = 0; |
| } |
| for (j = k + 1; j < n; j++) { |
| for (i = k + 1; i < m; i++) { |
| work[i] += e[j] * a[i][j]; |
| } |
| } |
| for (j = k + 1; j < n; j++) { |
| t = -e[j] / e[k + 1]; |
| for (i = k + 1; i < m; i++) { |
| a[i][j] += t * work[i]; |
| } |
| } |
| } |
| if (wantv) { |
| for (i = k + 1; i < n; i++) { |
| V[i][k] = e[i]; |
| } |
| } |
| } |
| } |
| |
| p = Math.min(n, m + 1); |
| if (nct < n) { |
| s[nct] = a[nct][nct]; |
| } |
| if (m < p) { |
| s[p - 1] = 0; |
| } |
| if (nrt + 1 < p) { |
| e[nrt] = a[nrt][p - 1]; |
| } |
| e[p - 1] = 0; |
| |
| if (wantu) { |
| for (j = nct; j < nu; j++) { |
| for (i = 0; i < m; i++) { |
| U[i][j] = 0; |
| } |
| U[j][j] = 1; |
| } |
| for (k = nct - 1; k >= 0; k--) { |
| if (s[k] !== 0) { |
| for (j = k + 1; j < nu; j++) { |
| t = 0; |
| for (i = k; i < m; i++) { |
| t += U[i][k] * U[i][j]; |
| } |
| t = -t / U[k][k]; |
| for (i = k; i < m; i++) { |
| U[i][j] += t * U[i][k]; |
| } |
| } |
| for (i = k; i < m; i++) { |
| U[i][k] = -U[i][k]; |
| } |
| U[k][k] = 1 + U[k][k]; |
| for (i = 0; i < k - 1; i++) { |
| U[i][k] = 0; |
| } |
| } else { |
| for (i = 0; i < m; i++) { |
| U[i][k] = 0; |
| } |
| U[k][k] = 1; |
| } |
| } |
| } |
| |
| if (wantv) { |
| for (k = n - 1; k >= 0; k--) { |
| if ((k < nrt) && (e[k] !== 0)) { |
| for (j = k + 1; j < n; j++) { |
| t = 0; |
| for (i = k + 1; i < n; i++) { |
| t += V[i][k] * V[i][j]; |
| } |
| t = -t / V[k + 1][k]; |
| for (i = k + 1; i < n; i++) { |
| V[i][j] += t * V[i][k]; |
| } |
| } |
| } |
| for (i = 0; i < n; i++) { |
| V[i][k] = 0; |
| } |
| V[k][k] = 1; |
| } |
| } |
| |
| var pp = p - 1, |
| iter = 0, |
| eps = Math.pow(2, -52); |
| while (p > 0) { |
| for (k = p - 2; k >= -1; k--) { |
| if (k === -1) { |
| break; |
| } |
| if (Math.abs(e[k]) <= eps * (Math.abs(s[k]) + Math.abs(s[k + 1]))) { |
| e[k] = 0; |
| break; |
| } |
| } |
| if (k === p - 2) { |
| kase = 4; |
| } else { |
| for (ks = p - 1; ks >= k; ks--) { |
| if (ks === k) { |
| break; |
| } |
| t = (ks !== p ? Math.abs(e[ks]) : 0) + (ks !== k + 1 ? Math.abs(e[ks - 1]) : 0); |
| if (Math.abs(s[ks]) <= eps * t) { |
| s[ks] = 0; |
| break; |
| } |
| } |
| if (ks === k) { |
| kase = 3; |
| } else if (ks === p - 1) { |
| kase = 1; |
| } else { |
| kase = 2; |
| k = ks; |
| } |
| } |
| |
| k++; |
| |
| switch (kase) { |
| case 1: { |
| f = e[p - 2]; |
| e[p - 2] = 0; |
| for (j = p - 2; j >= k; j--) { |
| t = hypotenuse(s[j], f); |
| cs = s[j] / t; |
| sn = f / t; |
| s[j] = t; |
| if (j !== k) { |
| f = -sn * e[j - 1]; |
| e[j - 1] = cs * e[j - 1]; |
| } |
| if (wantv) { |
| for (i = 0; i < n; i++) { |
| t = cs * V[i][j] + sn * V[i][p - 1]; |
| V[i][p - 1] = -sn * V[i][j] + cs * V[i][p - 1]; |
| V[i][j] = t; |
| } |
| } |
| } |
| break; |
| } |
| case 2 : { |
| f = e[k - 1]; |
| e[k - 1] = 0; |
| for (j = k; j < p; j++) { |
| t = hypotenuse(s[j], f); |
| cs = s[j] / t; |
| sn = f / t; |
| s[j] = t; |
| f = -sn * e[j]; |
| e[j] = cs * e[j]; |
| if (wantu) { |
| for (i = 0; i < m; i++) { |
| t = cs * U[i][j] + sn * U[i][k - 1]; |
| U[i][k - 1] = -sn * U[i][j] + cs * U[i][k - 1]; |
| U[i][j] = t; |
| } |
| } |
| } |
| break; |
| } |
| case 3 : { |
| scale = Math.max(Math.max(Math.max(Math.max(Math.abs(s[p - 1]), Math.abs(s[p - 2])), Math.abs(e[p - 2])), Math.abs(s[k])), Math.abs(e[k])); |
| sp = s[p - 1] / scale; |
| spm1 = s[p - 2] / scale; |
| epm1 = e[p - 2] / scale; |
| sk = s[k] / scale; |
| ek = e[k] / scale; |
| b = ((spm1 + sp) * (spm1 - sp) + epm1 * epm1) / 2; |
| c = (sp * epm1) * (sp * epm1); |
| shift = 0; |
| if ((b !== 0) || (c !== 0)) { |
| shift = Math.sqrt(b * b + c); |
| if (b < 0) { |
| shift = -shift; |
| } |
| shift = c / (b + shift); |
| } |
| f = (sk + sp) * (sk - sp) + shift; |
| g = sk * ek; |
| for (j = k; j < p - 1; j++) { |
| t = hypotenuse(f, g); |
| cs = f / t; |
| sn = g / t; |
| if (j !== k) { |
| e[j - 1] = t; |
| } |
| f = cs * s[j] + sn * e[j]; |
| e[j] = cs * e[j] - sn * s[j]; |
| g = sn * s[j + 1]; |
| s[j + 1] = cs * s[j + 1]; |
| if (wantv) { |
| for (i = 0; i < n; i++) { |
| t = cs * V[i][j] + sn * V[i][j + 1]; |
| V[i][j + 1] = -sn * V[i][j] + cs * V[i][j + 1]; |
| V[i][j] = t; |
| } |
| } |
| t = hypotenuse(f, g); |
| cs = f / t; |
| sn = g / t; |
| s[j] = t; |
| f = cs * e[j] + sn * s[j + 1]; |
| s[j + 1] = -sn * e[j] + cs * s[j + 1]; |
| g = sn * e[j + 1]; |
| e[j + 1] = cs * e[j + 1]; |
| if (wantu && (j < m - 1)) { |
| for (i = 0; i < m; i++) { |
| t = cs * U[i][j] + sn * U[i][j + 1]; |
| U[i][j + 1] = -sn * U[i][j] + cs * U[i][j + 1]; |
| U[i][j] = t; |
| } |
| } |
| } |
| e[p - 2] = f; |
| iter = iter + 1; |
| break; |
| } |
| case 4: { |
| if (s[k] <= 0) { |
| s[k] = (s[k] < 0 ? -s[k] : 0); |
| if (wantv) { |
| for (i = 0; i <= pp; i++) { |
| V[i][k] = -V[i][k]; |
| } |
| } |
| } |
| while (k < pp) { |
| if (s[k] >= s[k + 1]) { |
| break; |
| } |
| t = s[k]; |
| s[k] = s[k + 1]; |
| s[k + 1] = t; |
| if (wantv && (k < n - 1)) { |
| for (i = 0; i < n; i++) { |
| t = V[i][k + 1]; |
| V[i][k + 1] = V[i][k]; |
| V[i][k] = t; |
| } |
| } |
| if (wantu && (k < m - 1)) { |
| for (i = 0; i < m; i++) { |
| t = U[i][k + 1]; |
| U[i][k + 1] = U[i][k]; |
| U[i][k] = t; |
| } |
| } |
| k++; |
| } |
| iter = 0; |
| p--; |
| break; |
| } |
| // no default |
| } |
| } |
| |
| if (swapped) { |
| var tmp = V; |
| V = U; |
| U = tmp; |
| } |
| |
| this.m = m; |
| this.n = n; |
| this.s = s; |
| this.U = U; |
| this.V = V; |
| } |
| |
| SingularValueDecomposition.prototype = { |
| get condition() { |
| return this.s[0] / this.s[Math.min(this.m, this.n) - 1]; |
| }, |
| get norm2() { |
| return this.s[0]; |
| }, |
| get rank() { |
| var eps = Math.pow(2, -52), |
| tol = Math.max(this.m, this.n) * this.s[0] * eps, |
| r = 0, |
| s = this.s; |
| for (var i = 0, ii = s.length; i < ii; i++) { |
| if (s[i] > tol) { |
| r++; |
| } |
| } |
| return r; |
| }, |
| get diagonal() { |
| return this.s; |
| }, |
| // https://github.com/accord-net/framework/blob/development/Sources/Accord.Math/Decompositions/SingularValueDecomposition.cs |
| get threshold() { |
| return (Math.pow(2, -52) / 2) * Math.max(this.m, this.n) * this.s[0]; |
| }, |
| get leftSingularVectors() { |
| if (!Matrix.Matrix.isMatrix(this.U)) { |
| this.U = new Matrix.Matrix(this.U); |
| } |
| return this.U; |
| }, |
| get rightSingularVectors() { |
| if (!Matrix.Matrix.isMatrix(this.V)) { |
| this.V = new Matrix.Matrix(this.V); |
| } |
| return this.V; |
| }, |
| get diagonalMatrix() { |
| return Matrix.Matrix.diag(this.s); |
| }, |
| solve: function (value) { |
| |
| var Y = value, |
| e = this.threshold, |
| scols = this.s.length, |
| Ls = Matrix.Matrix.zeros(scols, scols), |
| i; |
| |
| for (i = 0; i < scols; i++) { |
| if (Math.abs(this.s[i]) <= e) { |
| Ls[i][i] = 0; |
| } else { |
| Ls[i][i] = 1 / this.s[i]; |
| } |
| } |
| |
| var U = this.U; |
| var V = this.rightSingularVectors; |
| |
| var VL = V.mmul(Ls), |
| vrows = V.rows, |
| urows = U.length, |
| VLU = Matrix.Matrix.zeros(vrows, urows), |
| j, k, sum; |
| |
| for (i = 0; i < vrows; i++) { |
| for (j = 0; j < urows; j++) { |
| sum = 0; |
| for (k = 0; k < scols; k++) { |
| sum += VL[i][k] * U[j][k]; |
| } |
| VLU[i][j] = sum; |
| } |
| } |
| |
| return VLU.mmul(Y); |
| }, |
| solveForDiagonal: function (value) { |
| return this.solve(Matrix.Matrix.diag(value)); |
| }, |
| inverse: function () { |
| var V = this.V; |
| var e = this.threshold, |
| vrows = V.length, |
| vcols = V[0].length, |
| X = new Matrix.Matrix(vrows, this.s.length), |
| i, j; |
| |
| for (i = 0; i < vrows; i++) { |
| for (j = 0; j < vcols; j++) { |
| if (Math.abs(this.s[j]) > e) { |
| X[i][j] = V[i][j] / this.s[j]; |
| } else { |
| X[i][j] = 0; |
| } |
| } |
| } |
| |
| var U = this.U; |
| |
| var urows = U.length, |
| ucols = U[0].length, |
| Y = new Matrix.Matrix(vrows, urows), |
| k, sum; |
| |
| for (i = 0; i < vrows; i++) { |
| for (j = 0; j < urows; j++) { |
| sum = 0; |
| for (k = 0; k < ucols; k++) { |
| sum += X[i][k] * U[j][k]; |
| } |
| Y[i][j] = sum; |
| } |
| } |
| |
| return Y; |
| } |
| }; |
| |
| MLMatrixDCSVD = SingularValueDecomposition; |
| } |
| |
| // ml-matrix src/abstractMatrix.js |
| let MLMatrixAbstractMatrix; |
| { |
| let LuDecomposition = MLMatrixDCLU; |
| let SvDecomposition = MLMatrixDCSVD; |
| let arrayUtils = MLArrayUtils; |
| let util = MLMatrixUtil; |
| |
| MLMatrixAbstractMatrix = function abstractMatrix(superCtor) { |
| if (superCtor === undefined) superCtor = Object; |
| |
| /** |
| * Real matrix |
| * @class Matrix |
| * @param {number|Array|Matrix} nRows - Number of rows of the new matrix, |
| * 2D array containing the data or Matrix instance to clone |
| * @param {number} [nColumns] - Number of columns of the new matrix |
| */ |
| class Matrix extends superCtor { |
| static get [Symbol.species]() { |
| return this; |
| } |
| |
| /** |
| * Constructs a Matrix with the chosen dimensions from a 1D array |
| * @param {number} newRows - Number of rows |
| * @param {number} newColumns - Number of columns |
| * @param {Array} newData - A 1D array containing data for the matrix |
| * @return {Matrix} - The new matrix |
| */ |
| static from1DArray(newRows, newColumns, newData) { |
| var length = newRows * newColumns; |
| if (length !== newData.length) { |
| throw new RangeError('Data length does not match given dimensions'); |
| } |
| var newMatrix = new this(newRows, newColumns); |
| for (var row = 0; row < newRows; row++) { |
| for (var column = 0; column < newColumns; column++) { |
| newMatrix.set(row, column, newData[row * newColumns + column]); |
| } |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Creates a row vector, a matrix with only one row. |
| * @param {Array} newData - A 1D array containing data for the vector |
| * @return {Matrix} - The new matrix |
| */ |
| static rowVector(newData) { |
| var vector = new this(1, newData.length); |
| for (var i = 0; i < newData.length; i++) { |
| vector.set(0, i, newData[i]); |
| } |
| return vector; |
| } |
| |
| /** |
| * Creates a column vector, a matrix with only one column. |
| * @param {Array} newData - A 1D array containing data for the vector |
| * @return {Matrix} - The new matrix |
| */ |
| static columnVector(newData) { |
| var vector = new this(newData.length, 1); |
| for (var i = 0; i < newData.length; i++) { |
| vector.set(i, 0, newData[i]); |
| } |
| return vector; |
| } |
| |
| /** |
| * Creates an empty matrix with the given dimensions. Values will be undefined. Same as using new Matrix(rows, columns). |
| * @param {number} rows - Number of rows |
| * @param {number} columns - Number of columns |
| * @return {Matrix} - The new matrix |
| */ |
| static empty(rows, columns) { |
| return new this(rows, columns); |
| } |
| |
| /** |
| * Creates a matrix with the given dimensions. Values will be set to zero. |
| * @param {number} rows - Number of rows |
| * @param {number} columns - Number of columns |
| * @return {Matrix} - The new matrix |
| */ |
| static zeros(rows, columns) { |
| return this.empty(rows, columns).fill(0); |
| } |
| |
| /** |
| * Creates a matrix with the given dimensions. Values will be set to one. |
| * @param {number} rows - Number of rows |
| * @param {number} columns - Number of columns |
| * @return {Matrix} - The new matrix |
| */ |
| static ones(rows, columns) { |
| return this.empty(rows, columns).fill(1); |
| } |
| |
| /** |
| * Creates a matrix with the given dimensions. Values will be randomly set. |
| * @param {number} rows - Number of rows |
| * @param {number} columns - Number of columns |
| * @param {function} [rng=Math.random] - Random number generator |
| * @return {Matrix} The new matrix |
| */ |
| static rand(rows, columns, rng) { |
| if (rng === undefined) rng = Math.random; |
| var matrix = this.empty(rows, columns); |
| for (var i = 0; i < rows; i++) { |
| for (var j = 0; j < columns; j++) { |
| matrix.set(i, j, rng()); |
| } |
| } |
| return matrix; |
| } |
| |
| /** |
| * Creates a matrix with the given dimensions. Values will be random integers. |
| * @param {number} rows - Number of rows |
| * @param {number} columns - Number of columns |
| * @param {number} [maxValue=1000] - Maximum value |
| * @param {function} [rng=Math.random] - Random number generator |
| * @return {Matrix} The new matrix |
| */ |
| static randInt(rows, columns, maxValue, rng) { |
| if (maxValue === undefined) maxValue = 1000; |
| if (rng === undefined) rng = Math.random; |
| var matrix = this.empty(rows, columns); |
| for (var i = 0; i < rows; i++) { |
| for (var j = 0; j < columns; j++) { |
| var value = Math.floor(rng() * maxValue); |
| matrix.set(i, j, value); |
| } |
| } |
| return matrix; |
| } |
| |
| /** |
| * Creates an identity matrix with the given dimension. Values of the diagonal will be 1 and others will be 0. |
| * @param {number} rows - Number of rows |
| * @param {number} [columns=rows] - Number of columns |
| * @param {number} [value=1] - Value to fill the diagonal with |
| * @return {Matrix} - The new identity matrix |
| */ |
| static eye(rows, columns, value) { |
| if (columns === undefined) columns = rows; |
| if (value === undefined) value = 1; |
| var min = Math.min(rows, columns); |
| var matrix = this.zeros(rows, columns); |
| for (var i = 0; i < min; i++) { |
| matrix.set(i, i, value); |
| } |
| return matrix; |
| } |
| |
| /** |
| * Creates a diagonal matrix based on the given array. |
| * @param {Array} data - Array containing the data for the diagonal |
| * @param {number} [rows] - Number of rows (Default: data.length) |
| * @param {number} [columns] - Number of columns (Default: rows) |
| * @return {Matrix} - The new diagonal matrix |
| */ |
| static diag(data, rows, columns) { |
| var l = data.length; |
| if (rows === undefined) rows = l; |
| if (columns === undefined) columns = rows; |
| var min = Math.min(l, rows, columns); |
| var matrix = this.zeros(rows, columns); |
| for (var i = 0; i < min; i++) { |
| matrix.set(i, i, data[i]); |
| } |
| return matrix; |
| } |
| |
| /** |
| * Returns a matrix whose elements are the minimum between matrix1 and matrix2 |
| * @param {Matrix} matrix1 |
| * @param {Matrix} matrix2 |
| * @return {Matrix} |
| */ |
| static min(matrix1, matrix2) { |
| matrix1 = this.checkMatrix(matrix1); |
| matrix2 = this.checkMatrix(matrix2); |
| var rows = matrix1.rows; |
| var columns = matrix1.columns; |
| var result = new this(rows, columns); |
| for (var i = 0; i < rows; i++) { |
| for (var j = 0; j < columns; j++) { |
| result.set(i, j, Math.min(matrix1.get(i, j), matrix2.get(i, j))); |
| } |
| } |
| return result; |
| } |
| |
| /** |
| * Returns a matrix whose elements are the maximum between matrix1 and matrix2 |
| * @param {Matrix} matrix1 |
| * @param {Matrix} matrix2 |
| * @return {Matrix} |
| */ |
| static max(matrix1, matrix2) { |
| matrix1 = this.checkMatrix(matrix1); |
| matrix2 = this.checkMatrix(matrix2); |
| var rows = matrix1.rows; |
| var columns = matrix1.columns; |
| var result = new this(rows, columns); |
| for (var i = 0; i < rows; i++) { |
| for (var j = 0; j < columns; j++) { |
| result.set(i, j, Math.max(matrix1.get(i, j), matrix2.get(i, j))); |
| } |
| } |
| return result; |
| } |
| |
| /** |
| * Check that the provided value is a Matrix and tries to instantiate one if not |
| * @param {*} value - The value to check |
| * @return {Matrix} |
| */ |
| static checkMatrix(value) { |
| return Matrix.isMatrix(value) ? value : new this(value); |
| } |
| |
| /** |
| * Returns true if the argument is a Matrix, false otherwise |
| * @param {*} value - The value to check |
| * @return {boolean} |
| */ |
| static isMatrix(value) { |
| return (value != null) && (value.klass === 'Matrix'); |
| } |
| |
| /** |
| * @prop {number} size - The number of elements in the matrix. |
| */ |
| get size() { |
| return this.rows * this.columns; |
| } |
| |
| /** |
| * Applies a callback for each element of the matrix. The function is called in the matrix (this) context. |
| * @param {function} callback - Function that will be called with two parameters : i (row) and j (column) |
| * @return {Matrix} this |
| */ |
| apply(callback) { |
| if (typeof callback !== 'function') { |
| throw new TypeError('callback must be a function'); |
| } |
| var ii = this.rows; |
| var jj = this.columns; |
| for (var i = 0; i < ii; i++) { |
| for (var j = 0; j < jj; j++) { |
| callback.call(this, i, j); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Returns a new 1D array filled row by row with the matrix values |
| * @return {Array} |
| */ |
| to1DArray() { |
| var array = new Array(this.size); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| array[i * this.columns + j] = this.get(i, j); |
| } |
| } |
| return array; |
| } |
| |
| /** |
| * Returns a 2D array containing a copy of the data |
| * @return {Array} |
| */ |
| to2DArray() { |
| var copy = new Array(this.rows); |
| for (var i = 0; i < this.rows; i++) { |
| copy[i] = new Array(this.columns); |
| for (var j = 0; j < this.columns; j++) { |
| copy[i][j] = this.get(i, j); |
| } |
| } |
| return copy; |
| } |
| |
| /** |
| * @return {boolean} true if the matrix has one row |
| */ |
| isRowVector() { |
| return this.rows === 1; |
| } |
| |
| /** |
| * @return {boolean} true if the matrix has one column |
| */ |
| isColumnVector() { |
| return this.columns === 1; |
| } |
| |
| /** |
| * @return {boolean} true if the matrix has one row or one column |
| */ |
| isVector() { |
| return (this.rows === 1) || (this.columns === 1); |
| } |
| |
| /** |
| * @return {boolean} true if the matrix has the same number of rows and columns |
| */ |
| isSquare() { |
| return this.rows === this.columns; |
| } |
| |
| /** |
| * @return {boolean} true if the matrix is square and has the same values on both sides of the diagonal |
| */ |
| isSymmetric() { |
| if (this.isSquare()) { |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j <= i; j++) { |
| if (this.get(i, j) !== this.get(j, i)) { |
| return false; |
| } |
| } |
| } |
| return true; |
| } |
| return false; |
| } |
| |
| /** |
| * Sets a given element of the matrix. mat.set(3,4,1) is equivalent to mat[3][4]=1 |
| * @abstract |
| * @param {number} rowIndex - Index of the row |
| * @param {number} columnIndex - Index of the column |
| * @param {number} value - The new value for the element |
| * @return {Matrix} this |
| */ |
| set(rowIndex, columnIndex, value) { // eslint-disable-line no-unused-vars |
| throw new Error('set method is unimplemented'); |
| } |
| |
| /** |
| * Returns the given element of the matrix. mat.get(3,4) is equivalent to matrix[3][4] |
| * @abstract |
| * @param {number} rowIndex - Index of the row |
| * @param {number} columnIndex - Index of the column |
| * @return {number} |
| */ |
| get(rowIndex, columnIndex) { // eslint-disable-line no-unused-vars |
| throw new Error('get method is unimplemented'); |
| } |
| |
| /** |
| * Creates a new matrix that is a repetition of the current matrix. New matrix has rowRep times the number of |
| * rows of the matrix, and colRep times the number of columns of the matrix |
| * @param {number} rowRep - Number of times the rows should be repeated |
| * @param {number} colRep - Number of times the columns should be re |
| * @return {Matrix} |
| * @example |
| * var matrix = new Matrix([[1,2]]); |
| * matrix.repeat(2); // [[1,2],[1,2]] |
| */ |
| repeat(rowRep, colRep) { |
| rowRep = rowRep || 1; |
| colRep = colRep || 1; |
| var matrix = new this.constructor[Symbol.species](this.rows * rowRep, this.columns * colRep); |
| for (var i = 0; i < rowRep; i++) { |
| for (var j = 0; j < colRep; j++) { |
| matrix.setSubMatrix(this, this.rows * i, this.columns * j); |
| } |
| } |
| return matrix; |
| } |
| |
| /** |
| * Fills the matrix with a given value. All elements will be set to this value. |
| * @param {number} value - New value |
| * @return {Matrix} this |
| */ |
| fill(value) { |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, value); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Negates the matrix. All elements will be multiplied by (-1) |
| * @return {Matrix} this |
| */ |
| neg() { |
| return this.mulS(-1); |
| } |
| |
| /** |
| * Returns a new array from the given row index |
| * @param {number} index - Row index |
| * @return {Array} |
| */ |
| getRow(index) { |
| util.checkRowIndex(this, index); |
| var row = new Array(this.columns); |
| for (var i = 0; i < this.columns; i++) { |
| row[i] = this.get(index, i); |
| } |
| return row; |
| } |
| |
| /** |
| * Returns a new row vector from the given row index |
| * @param {number} index - Row index |
| * @return {Matrix} |
| */ |
| getRowVector(index) { |
| return this.constructor.rowVector(this.getRow(index)); |
| } |
| |
| /** |
| * Sets a row at the given index |
| * @param {number} index - Row index |
| * @param {Array|Matrix} array - Array or vector |
| * @return {Matrix} this |
| */ |
| setRow(index, array) { |
| util.checkRowIndex(this, index); |
| array = util.checkRowVector(this, array); |
| for (var i = 0; i < this.columns; i++) { |
| this.set(index, i, array[i]); |
| } |
| return this; |
| } |
| |
| /** |
| * Swaps two rows |
| * @param {number} row1 - First row index |
| * @param {number} row2 - Second row index |
| * @return {Matrix} this |
| */ |
| swapRows(row1, row2) { |
| util.checkRowIndex(this, row1); |
| util.checkRowIndex(this, row2); |
| for (var i = 0; i < this.columns; i++) { |
| var temp = this.get(row1, i); |
| this.set(row1, i, this.get(row2, i)); |
| this.set(row2, i, temp); |
| } |
| return this; |
| } |
| |
| /** |
| * Returns a new array from the given column index |
| * @param {number} index - Column index |
| * @return {Array} |
| */ |
| getColumn(index) { |
| util.checkColumnIndex(this, index); |
| var column = new Array(this.rows); |
| for (var i = 0; i < this.rows; i++) { |
| column[i] = this.get(i, index); |
| } |
| return column; |
| } |
| |
| /** |
| * Returns a new column vector from the given column index |
| * @param {number} index - Column index |
| * @return {Matrix} |
| */ |
| getColumnVector(index) { |
| return this.constructor.columnVector(this.getColumn(index)); |
| } |
| |
| /** |
| * Sets a column at the given index |
| * @param {number} index - Column index |
| * @param {Array|Matrix} array - Array or vector |
| * @return {Matrix} this |
| */ |
| setColumn(index, array) { |
| util.checkColumnIndex(this, index); |
| array = util.checkColumnVector(this, array); |
| for (var i = 0; i < this.rows; i++) { |
| this.set(i, index, array[i]); |
| } |
| return this; |
| } |
| |
| /** |
| * Swaps two columns |
| * @param {number} column1 - First column index |
| * @param {number} column2 - Second column index |
| * @return {Matrix} this |
| */ |
| swapColumns(column1, column2) { |
| util.checkColumnIndex(this, column1); |
| util.checkColumnIndex(this, column2); |
| for (var i = 0; i < this.rows; i++) { |
| var temp = this.get(i, column1); |
| this.set(i, column1, this.get(i, column2)); |
| this.set(i, column2, temp); |
| } |
| return this; |
| } |
| |
| /** |
| * Adds the values of a vector to each row |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| addRowVector(vector) { |
| vector = util.checkRowVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) + vector[j]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Subtracts the values of a vector from each row |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| subRowVector(vector) { |
| vector = util.checkRowVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) - vector[j]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Multiplies the values of a vector with each row |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| mulRowVector(vector) { |
| vector = util.checkRowVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) * vector[j]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Divides the values of each row by those of a vector |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| divRowVector(vector) { |
| vector = util.checkRowVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) / vector[j]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Adds the values of a vector to each column |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| addColumnVector(vector) { |
| vector = util.checkColumnVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) + vector[i]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Subtracts the values of a vector from each column |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| subColumnVector(vector) { |
| vector = util.checkColumnVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) - vector[i]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Multiplies the values of a vector with each column |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| mulColumnVector(vector) { |
| vector = util.checkColumnVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) * vector[i]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Divides the values of each column by those of a vector |
| * @param {Array|Matrix} vector - Array or vector |
| * @return {Matrix} this |
| */ |
| divColumnVector(vector) { |
| vector = util.checkColumnVector(this, vector); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) / vector[i]); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Multiplies the values of a row with a scalar |
| * @param {number} index - Row index |
| * @param {number} value |
| * @return {Matrix} this |
| */ |
| mulRow(index, value) { |
| util.checkRowIndex(this, index); |
| for (var i = 0; i < this.columns; i++) { |
| this.set(index, i, this.get(index, i) * value); |
| } |
| return this; |
| } |
| |
| /** |
| * Multiplies the values of a column with a scalar |
| * @param {number} index - Column index |
| * @param {number} value |
| * @return {Matrix} this |
| */ |
| mulColumn(index, value) { |
| util.checkColumnIndex(this, index); |
| for (var i = 0; i < this.rows; i++) { |
| this.set(i, index, this.get(i, index) * value); |
| } |
| return this; |
| } |
| |
| /** |
| * Returns the maximum value of the matrix |
| * @return {number} |
| */ |
| max() { |
| var v = this.get(0, 0); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| if (this.get(i, j) > v) { |
| v = this.get(i, j); |
| } |
| } |
| } |
| return v; |
| } |
| |
| /** |
| * Returns the index of the maximum value |
| * @return {Array} |
| */ |
| maxIndex() { |
| var v = this.get(0, 0); |
| var idx = [0, 0]; |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| if (this.get(i, j) > v) { |
| v = this.get(i, j); |
| idx[0] = i; |
| idx[1] = j; |
| } |
| } |
| } |
| return idx; |
| } |
| |
| /** |
| * Returns the minimum value of the matrix |
| * @return {number} |
| */ |
| min() { |
| var v = this.get(0, 0); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| if (this.get(i, j) < v) { |
| v = this.get(i, j); |
| } |
| } |
| } |
| return v; |
| } |
| |
| /** |
| * Returns the index of the minimum value |
| * @return {Array} |
| */ |
| minIndex() { |
| var v = this.get(0, 0); |
| var idx = [0, 0]; |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| if (this.get(i, j) < v) { |
| v = this.get(i, j); |
| idx[0] = i; |
| idx[1] = j; |
| } |
| } |
| } |
| return idx; |
| } |
| |
| /** |
| * Returns the maximum value of one row |
| * @param {number} row - Row index |
| * @return {number} |
| */ |
| maxRow(row) { |
| util.checkRowIndex(this, row); |
| var v = this.get(row, 0); |
| for (var i = 1; i < this.columns; i++) { |
| if (this.get(row, i) > v) { |
| v = this.get(row, i); |
| } |
| } |
| return v; |
| } |
| |
| /** |
| * Returns the index of the maximum value of one row |
| * @param {number} row - Row index |
| * @return {Array} |
| */ |
| maxRowIndex(row) { |
| util.checkRowIndex(this, row); |
| var v = this.get(row, 0); |
| var idx = [row, 0]; |
| for (var i = 1; i < this.columns; i++) { |
| if (this.get(row, i) > v) { |
| v = this.get(row, i); |
| idx[1] = i; |
| } |
| } |
| return idx; |
| } |
| |
| /** |
| * Returns the minimum value of one row |
| * @param {number} row - Row index |
| * @return {number} |
| */ |
| minRow(row) { |
| util.checkRowIndex(this, row); |
| var v = this.get(row, 0); |
| for (var i = 1; i < this.columns; i++) { |
| if (this.get(row, i) < v) { |
| v = this.get(row, i); |
| } |
| } |
| return v; |
| } |
| |
| /** |
| * Returns the index of the maximum value of one row |
| * @param {number} row - Row index |
| * @return {Array} |
| */ |
| minRowIndex(row) { |
| util.checkRowIndex(this, row); |
| var v = this.get(row, 0); |
| var idx = [row, 0]; |
| for (var i = 1; i < this.columns; i++) { |
| if (this.get(row, i) < v) { |
| v = this.get(row, i); |
| idx[1] = i; |
| } |
| } |
| return idx; |
| } |
| |
| /** |
| * Returns the maximum value of one column |
| * @param {number} column - Column index |
| * @return {number} |
| */ |
| maxColumn(column) { |
| util.checkColumnIndex(this, column); |
| var v = this.get(0, column); |
| for (var i = 1; i < this.rows; i++) { |
| if (this.get(i, column) > v) { |
| v = this.get(i, column); |
| } |
| } |
| return v; |
| } |
| |
| /** |
| * Returns the index of the maximum value of one column |
| * @param {number} column - Column index |
| * @return {Array} |
| */ |
| maxColumnIndex(column) { |
| util.checkColumnIndex(this, column); |
| var v = this.get(0, column); |
| var idx = [0, column]; |
| for (var i = 1; i < this.rows; i++) { |
| if (this.get(i, column) > v) { |
| v = this.get(i, column); |
| idx[0] = i; |
| } |
| } |
| return idx; |
| } |
| |
| /** |
| * Returns the minimum value of one column |
| * @param {number} column - Column index |
| * @return {number} |
| */ |
| minColumn(column) { |
| util.checkColumnIndex(this, column); |
| var v = this.get(0, column); |
| for (var i = 1; i < this.rows; i++) { |
| if (this.get(i, column) < v) { |
| v = this.get(i, column); |
| } |
| } |
| return v; |
| } |
| |
| /** |
| * Returns the index of the minimum value of one column |
| * @param {number} column - Column index |
| * @return {Array} |
| */ |
| minColumnIndex(column) { |
| util.checkColumnIndex(this, column); |
| var v = this.get(0, column); |
| var idx = [0, column]; |
| for (var i = 1; i < this.rows; i++) { |
| if (this.get(i, column) < v) { |
| v = this.get(i, column); |
| idx[0] = i; |
| } |
| } |
| return idx; |
| } |
| |
| /** |
| * Returns an array containing the diagonal values of the matrix |
| * @return {Array} |
| */ |
| diag() { |
| var min = Math.min(this.rows, this.columns); |
| var diag = new Array(min); |
| for (var i = 0; i < min; i++) { |
| diag[i] = this.get(i, i); |
| } |
| return diag; |
| } |
| |
| /** |
| * Returns the sum by the argument given, if no argument given, |
| * it returns the sum of all elements of the matrix. |
| * @param {string} by - sum by 'row' or 'column'. |
| * @return {Matrix|number} |
| */ |
| sum(by) { |
| switch (by) { |
| case 'row': |
| return util.sumByRow(this); |
| case 'column': |
| return util.sumByColumn(this); |
| default: |
| return util.sumAll(this); |
| } |
| } |
| |
| /** |
| * Returns the mean of all elements of the matrix |
| * @return {number} |
| */ |
| mean() { |
| return this.sum() / this.size; |
| } |
| |
| /** |
| * Returns the product of all elements of the matrix |
| * @return {number} |
| */ |
| prod() { |
| var prod = 1; |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| prod *= this.get(i, j); |
| } |
| } |
| return prod; |
| } |
| |
| /** |
| * Computes the cumulative sum of the matrix elements (in place, row by row) |
| * @return {Matrix} this |
| */ |
| cumulativeSum() { |
| var sum = 0; |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| sum += this.get(i, j); |
| this.set(i, j, sum); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Computes the dot (scalar) product between the matrix and another |
| * @param {Matrix} vector2 vector |
| * @return {number} |
| */ |
| dot(vector2) { |
| if (Matrix.isMatrix(vector2)) vector2 = vector2.to1DArray(); |
| var vector1 = this.to1DArray(); |
| if (vector1.length !== vector2.length) { |
| throw new RangeError('vectors do not have the same size'); |
| } |
| var dot = 0; |
| for (var i = 0; i < vector1.length; i++) { |
| dot += vector1[i] * vector2[i]; |
| } |
| return dot; |
| } |
| |
| /** |
| * Returns the matrix product between this and other |
| * @param {Matrix} other |
| * @return {Matrix} |
| */ |
| mmul(other) { |
| other = this.constructor.checkMatrix(other); |
| if (this.columns !== other.rows) { |
| // eslint-disable-next-line no-console |
| console.warn('Number of columns of left matrix are not equal to number of rows of right matrix.'); |
| } |
| |
| var m = this.rows; |
| var n = this.columns; |
| var p = other.columns; |
| |
| var result = new this.constructor[Symbol.species](m, p); |
| |
| var Bcolj = new Array(n); |
| for (var j = 0; j < p; j++) { |
| for (var k = 0; k < n; k++) { |
| Bcolj[k] = other.get(k, j); |
| } |
| |
| for (var i = 0; i < m; i++) { |
| var s = 0; |
| for (k = 0; k < n; k++) { |
| s += this.get(i, k) * Bcolj[k]; |
| } |
| |
| result.set(i, j, s); |
| } |
| } |
| return result; |
| } |
| |
| strassen2x2(other) { |
| var result = new this.constructor[Symbol.species](2, 2); |
| const a11 = this.get(0, 0); |
| const b11 = other.get(0, 0); |
| const a12 = this.get(0, 1); |
| const b12 = other.get(0, 1); |
| const a21 = this.get(1, 0); |
| const b21 = other.get(1, 0); |
| const a22 = this.get(1, 1); |
| const b22 = other.get(1, 1); |
| |
| // Compute intermediate values. |
| const m1 = (a11 + a22) * (b11 + b22); |
| const m2 = (a21 + a22) * b11; |
| const m3 = a11 * (b12 - b22); |
| const m4 = a22 * (b21 - b11); |
| const m5 = (a11 + a12) * b22; |
| const m6 = (a21 - a11) * (b11 + b12); |
| const m7 = (a12 - a22) * (b21 + b22); |
| |
| // Combine intermediate values into the output. |
| const c00 = m1 + m4 - m5 + m7; |
| const c01 = m3 + m5; |
| const c10 = m2 + m4; |
| const c11 = m1 - m2 + m3 + m6; |
| |
| result.set(0, 0, c00); |
| result.set(0, 1, c01); |
| result.set(1, 0, c10); |
| result.set(1, 1, c11); |
| return result; |
| } |
| |
| strassen3x3(other) { |
| var result = new this.constructor[Symbol.species](3, 3); |
| |
| const a00 = this.get(0, 0); |
| const a01 = this.get(0, 1); |
| const a02 = this.get(0, 2); |
| const a10 = this.get(1, 0); |
| const a11 = this.get(1, 1); |
| const a12 = this.get(1, 2); |
| const a20 = this.get(2, 0); |
| const a21 = this.get(2, 1); |
| const a22 = this.get(2, 2); |
| |
| const b00 = other.get(0, 0); |
| const b01 = other.get(0, 1); |
| const b02 = other.get(0, 2); |
| const b10 = other.get(1, 0); |
| const b11 = other.get(1, 1); |
| const b12 = other.get(1, 2); |
| const b20 = other.get(2, 0); |
| const b21 = other.get(2, 1); |
| const b22 = other.get(2, 2); |
| |
| const m1 = (a00 + a01 + a02 - a10 - a11 - a21 - a22) * b11; |
| const m2 = (a00 - a10) * (-b01 + b11); |
| const m3 = a11 * (-b00 + b01 + b10 - b11 - b12 - b20 + b22); |
| const m4 = (-a00 + a10 + a11) * (b00 - b01 + b11); |
| const m5 = (a10 + a11) * (-b00 + b01); |
| const m6 = a00 * b00; |
| const m7 = (-a00 + a20 + a21) * (b00 - b02 + b12); |
| const m8 = (-a00 + a20) * (b02 - b12); |
| const m9 = (a20 + a21) * (-b00 + b02); |
| const m10 = (a00 + a01 + a02 - a11 - a12 - a20 - a21) * b12; |
| const m11 = a21 * (-b00 + b02 + b10 - b11 - b12 - b20 + b21); |
| const m12 = (-a02 + a21 + a22) * (b11 + b20 - b21); |
| const m13 = (a02 - a22) * (b11 - b21); |
| const m14 = a02 * b20; |
| const m15 = (a21 + a22) * (-b20 + b21); |
| const m16 = (-a02 + a11 + a12) * (b12 + b20 - b22); |
| const m17 = (a02 - a12) * (b12 - b22); |
| const m18 = (a11 + a12) * (-b20 + b22); |
| const m19 = a01 * b10; |
| const m20 = a12 * b21; |
| const m21 = a10 * b02; |
| const m22 = a20 * b01; |
| const m23 = a22 * b22; |
| |
| const c00 = m6 + m14 + m19; |
| const c01 = m1 + m4 + m5 + m6 + m12 + m14 + m15; |
| const c02 = m6 + m7 + m9 + m10 + m14 + m16 + m18; |
| const c10 = m2 + m3 + m4 + m6 + m14 + m16 + m17; |
| const c11 = m2 + m4 + m5 + m6 + m20; |
| const c12 = m14 + m16 + m17 + m18 + m21; |
| const c20 = m6 + m7 + m8 + m11 + m12 + m13 + m14; |
| const c21 = m12 + m13 + m14 + m15 + m22; |
| const c22 = m6 + m7 + m8 + m9 + m23; |
| |
| result.set(0, 0, c00); |
| result.set(0, 1, c01); |
| result.set(0, 2, c02); |
| result.set(1, 0, c10); |
| result.set(1, 1, c11); |
| result.set(1, 2, c12); |
| result.set(2, 0, c20); |
| result.set(2, 1, c21); |
| result.set(2, 2, c22); |
| return result; |
| } |
| |
| /** |
| * Returns the matrix product between x and y. More efficient than mmul(other) only when we multiply squared matrix and when the size of the matrix is > 1000. |
| * @param {Matrix} y |
| * @return {Matrix} |
| */ |
| mmulStrassen(y) { |
| var x = this.clone(); |
| var r1 = x.rows; |
| var c1 = x.columns; |
| var r2 = y.rows; |
| var c2 = y.columns; |
| if (c1 !== r2) { |
| // eslint-disable-next-line no-console |
| console.warn(`Multiplying ${r1} x ${c1} and ${r2} x ${c2} matrix: dimensions do not match.`); |
| } |
| |
| // Put a matrix into the top left of a matrix of zeros. |
| // `rows` and `cols` are the dimensions of the output matrix. |
| function embed(mat, rows, cols) { |
| var r = mat.rows; |
| var c = mat.columns; |
| if ((r === rows) && (c === cols)) { |
| return mat; |
| } else { |
| var resultat = Matrix.zeros(rows, cols); |
| resultat = resultat.setSubMatrix(mat, 0, 0); |
| return resultat; |
| } |
| } |
| |
| |
| // Make sure both matrices are the same size. |
| // This is exclusively for simplicity: |
| // this algorithm can be implemented with matrices of different sizes. |
| |
| var r = Math.max(r1, r2); |
| var c = Math.max(c1, c2); |
| x = embed(x, r, c); |
| y = embed(y, r, c); |
| |
| // Our recursive multiplication function. |
| function blockMult(a, b, rows, cols) { |
| // For small matrices, resort to naive multiplication. |
| if (rows <= 512 || cols <= 512) { |
| return a.mmul(b); // a is equivalent to this |
| } |
| |
| // Apply dynamic padding. |
| if ((rows % 2 === 1) && (cols % 2 === 1)) { |
| a = embed(a, rows + 1, cols + 1); |
| b = embed(b, rows + 1, cols + 1); |
| } else if (rows % 2 === 1) { |
| a = embed(a, rows + 1, cols); |
| b = embed(b, rows + 1, cols); |
| } else if (cols % 2 === 1) { |
| a = embed(a, rows, cols + 1); |
| b = embed(b, rows, cols + 1); |
| } |
| |
| var halfRows = parseInt(a.rows / 2); |
| var halfCols = parseInt(a.columns / 2); |
| // Subdivide input matrices. |
| var a11 = a.subMatrix(0, halfRows - 1, 0, halfCols - 1); |
| var b11 = b.subMatrix(0, halfRows - 1, 0, halfCols - 1); |
| |
| var a12 = a.subMatrix(0, halfRows - 1, halfCols, a.columns - 1); |
| var b12 = b.subMatrix(0, halfRows - 1, halfCols, b.columns - 1); |
| |
| var a21 = a.subMatrix(halfRows, a.rows - 1, 0, halfCols - 1); |
| var b21 = b.subMatrix(halfRows, b.rows - 1, 0, halfCols - 1); |
| |
| var a22 = a.subMatrix(halfRows, a.rows - 1, halfCols, a.columns - 1); |
| var b22 = b.subMatrix(halfRows, b.rows - 1, halfCols, b.columns - 1); |
| |
| // Compute intermediate values. |
| var m1 = blockMult(Matrix.add(a11, a22), Matrix.add(b11, b22), halfRows, halfCols); |
| var m2 = blockMult(Matrix.add(a21, a22), b11, halfRows, halfCols); |
| var m3 = blockMult(a11, Matrix.sub(b12, b22), halfRows, halfCols); |
| var m4 = blockMult(a22, Matrix.sub(b21, b11), halfRows, halfCols); |
| var m5 = blockMult(Matrix.add(a11, a12), b22, halfRows, halfCols); |
| var m6 = blockMult(Matrix.sub(a21, a11), Matrix.add(b11, b12), halfRows, halfCols); |
| var m7 = blockMult(Matrix.sub(a12, a22), Matrix.add(b21, b22), halfRows, halfCols); |
| |
| // Combine intermediate values into the output. |
| var c11 = Matrix.add(m1, m4); |
| c11.sub(m5); |
| c11.add(m7); |
| var c12 = Matrix.add(m3, m5); |
| var c21 = Matrix.add(m2, m4); |
| var c22 = Matrix.sub(m1, m2); |
| c22.add(m3); |
| c22.add(m6); |
| |
| //Crop output to the desired size (undo dynamic padding). |
| var resultat = Matrix.zeros(2 * c11.rows, 2 * c11.columns); |
| resultat = resultat.setSubMatrix(c11, 0, 0); |
| resultat = resultat.setSubMatrix(c12, c11.rows, 0); |
| resultat = resultat.setSubMatrix(c21, 0, c11.columns); |
| resultat = resultat.setSubMatrix(c22, c11.rows, c11.columns); |
| return resultat.subMatrix(0, rows - 1, 0, cols - 1); |
| } |
| return blockMult(x, y, r, c); |
| } |
| |
| /** |
| * Returns a row-by-row scaled matrix |
| * @param {number} [min=0] - Minimum scaled value |
| * @param {number} [max=1] - Maximum scaled value |
| * @return {Matrix} - The scaled matrix |
| */ |
| scaleRows(min, max) { |
| min = min === undefined ? 0 : min; |
| max = max === undefined ? 1 : max; |
| if (min >= max) { |
| throw new RangeError('min should be strictly smaller than max'); |
| } |
| var newMatrix = this.constructor.empty(this.rows, this.columns); |
| for (var i = 0; i < this.rows; i++) { |
| var scaled = arrayUtils.scale(this.getRow(i), {min, max}); |
| newMatrix.setRow(i, scaled); |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Returns a new column-by-column scaled matrix |
| * @param {number} [min=0] - Minimum scaled value |
| * @param {number} [max=1] - Maximum scaled value |
| * @return {Matrix} - The new scaled matrix |
| * @example |
| * var matrix = new Matrix([[1,2],[-1,0]]); |
| * var scaledMatrix = matrix.scaleColumns(); // [[1,1],[0,0]] |
| */ |
| scaleColumns(min, max) { |
| min = min === undefined ? 0 : min; |
| max = max === undefined ? 1 : max; |
| if (min >= max) { |
| throw new RangeError('min should be strictly smaller than max'); |
| } |
| var newMatrix = this.constructor.empty(this.rows, this.columns); |
| for (var i = 0; i < this.columns; i++) { |
| var scaled = arrayUtils.scale(this.getColumn(i), { |
| min: min, |
| max: max |
| }); |
| newMatrix.setColumn(i, scaled); |
| } |
| return newMatrix; |
| } |
| |
| |
| /** |
| * Returns the Kronecker product (also known as tensor product) between this and other |
| * See https://en.wikipedia.org/wiki/Kronecker_product |
| * @param {Matrix} other |
| * @return {Matrix} |
| */ |
| kroneckerProduct(other) { |
| other = this.constructor.checkMatrix(other); |
| |
| var m = this.rows; |
| var n = this.columns; |
| var p = other.rows; |
| var q = other.columns; |
| |
| var result = new this.constructor[Symbol.species](m * p, n * q); |
| for (var i = 0; i < m; i++) { |
| for (var j = 0; j < n; j++) { |
| for (var k = 0; k < p; k++) { |
| for (var l = 0; l < q; l++) { |
| result[p * i + k][q * j + l] = this.get(i, j) * other.get(k, l); |
| } |
| } |
| } |
| } |
| return result; |
| } |
| |
| /** |
| * Transposes the matrix and returns a new one containing the result |
| * @return {Matrix} |
| */ |
| transpose() { |
| var result = new this.constructor[Symbol.species](this.columns, this.rows); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| result.set(j, i, this.get(i, j)); |
| } |
| } |
| return result; |
| } |
| |
| /** |
| * Sorts the rows (in place) |
| * @param {function} compareFunction - usual Array.prototype.sort comparison function |
| * @return {Matrix} this |
| */ |
| sortRows(compareFunction) { |
| if (compareFunction === undefined) compareFunction = compareNumbers; |
| for (var i = 0; i < this.rows; i++) { |
| this.setRow(i, this.getRow(i).sort(compareFunction)); |
| } |
| return this; |
| } |
| |
| /** |
| * Sorts the columns (in place) |
| * @param {function} compareFunction - usual Array.prototype.sort comparison function |
| * @return {Matrix} this |
| */ |
| sortColumns(compareFunction) { |
| if (compareFunction === undefined) compareFunction = compareNumbers; |
| for (var i = 0; i < this.columns; i++) { |
| this.setColumn(i, this.getColumn(i).sort(compareFunction)); |
| } |
| return this; |
| } |
| |
| /** |
| * Returns a subset of the matrix |
| * @param {number} startRow - First row index |
| * @param {number} endRow - Last row index |
| * @param {number} startColumn - First column index |
| * @param {number} endColumn - Last column index |
| * @return {Matrix} |
| */ |
| subMatrix(startRow, endRow, startColumn, endColumn) { |
| util.checkRange(this, startRow, endRow, startColumn, endColumn); |
| var newMatrix = new this.constructor[Symbol.species](endRow - startRow + 1, endColumn - startColumn + 1); |
| for (var i = startRow; i <= endRow; i++) { |
| for (var j = startColumn; j <= endColumn; j++) { |
| newMatrix[i - startRow][j - startColumn] = this.get(i, j); |
| } |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Returns a subset of the matrix based on an array of row indices |
| * @param {Array} indices - Array containing the row indices |
| * @param {number} [startColumn = 0] - First column index |
| * @param {number} [endColumn = this.columns-1] - Last column index |
| * @return {Matrix} |
| */ |
| subMatrixRow(indices, startColumn, endColumn) { |
| if (startColumn === undefined) startColumn = 0; |
| if (endColumn === undefined) endColumn = this.columns - 1; |
| if ((startColumn > endColumn) || (startColumn < 0) || (startColumn >= this.columns) || (endColumn < 0) || (endColumn >= this.columns)) { |
| throw new RangeError('Argument out of range'); |
| } |
| |
| var newMatrix = new this.constructor[Symbol.species](indices.length, endColumn - startColumn + 1); |
| for (var i = 0; i < indices.length; i++) { |
| for (var j = startColumn; j <= endColumn; j++) { |
| if (indices[i] < 0 || indices[i] >= this.rows) { |
| throw new RangeError('Row index out of range: ' + indices[i]); |
| } |
| newMatrix.set(i, j - startColumn, this.get(indices[i], j)); |
| } |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Returns a subset of the matrix based on an array of column indices |
| * @param {Array} indices - Array containing the column indices |
| * @param {number} [startRow = 0] - First row index |
| * @param {number} [endRow = this.rows-1] - Last row index |
| * @return {Matrix} |
| */ |
| subMatrixColumn(indices, startRow, endRow) { |
| if (startRow === undefined) startRow = 0; |
| if (endRow === undefined) endRow = this.rows - 1; |
| if ((startRow > endRow) || (startRow < 0) || (startRow >= this.rows) || (endRow < 0) || (endRow >= this.rows)) { |
| throw new RangeError('Argument out of range'); |
| } |
| |
| var newMatrix = new this.constructor[Symbol.species](endRow - startRow + 1, indices.length); |
| for (var i = 0; i < indices.length; i++) { |
| for (var j = startRow; j <= endRow; j++) { |
| if (indices[i] < 0 || indices[i] >= this.columns) { |
| throw new RangeError('Column index out of range: ' + indices[i]); |
| } |
| newMatrix.set(j - startRow, i, this.get(j, indices[i])); |
| } |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Set a part of the matrix to the given sub-matrix |
| * @param {Matrix|Array< Array >} matrix - The source matrix from which to extract values. |
| * @param {number} startRow - The index of the first row to set |
| * @param {number} startColumn - The index of the first column to set |
| * @return {Matrix} |
| */ |
| setSubMatrix(matrix, startRow, startColumn) { |
| matrix = this.constructor.checkMatrix(matrix); |
| var endRow = startRow + matrix.rows - 1; |
| var endColumn = startColumn + matrix.columns - 1; |
| util.checkRange(this, startRow, endRow, startColumn, endColumn); |
| for (var i = 0; i < matrix.rows; i++) { |
| for (var j = 0; j < matrix.columns; j++) { |
| this[startRow + i][startColumn + j] = matrix.get(i, j); |
| } |
| } |
| return this; |
| } |
| |
| /** |
| * Return a new matrix based on a selection of rows and columns |
| * @param {Array<number>} rowIndices - The row indices to select. Order matters and an index can be more than once. |
| * @param {Array<number>} columnIndices - The column indices to select. Order matters and an index can be use more than once. |
| * @return {Matrix} The new matrix |
| */ |
| selection(rowIndices, columnIndices) { |
| var indices = util.checkIndices(this, rowIndices, columnIndices); |
| var newMatrix = new this.constructor[Symbol.species](rowIndices.length, columnIndices.length); |
| for (var i = 0; i < indices.row.length; i++) { |
| var rowIndex = indices.row[i]; |
| for (var j = 0; j < indices.column.length; j++) { |
| var columnIndex = indices.column[j]; |
| newMatrix[i][j] = this.get(rowIndex, columnIndex); |
| } |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Returns the trace of the matrix (sum of the diagonal elements) |
| * @return {number} |
| */ |
| trace() { |
| var min = Math.min(this.rows, this.columns); |
| var trace = 0; |
| for (var i = 0; i < min; i++) { |
| trace += this.get(i, i); |
| } |
| return trace; |
| } |
| |
| /* |
| Matrix views |
| */ |
| |
| /** |
| * Returns a view of the transposition of the matrix |
| * @return {MatrixTransposeView} |
| */ |
| transposeView() { |
| return new MLMatrixTransposeView(this); |
| } |
| |
| /** |
| * Returns a view of the row vector with the given index |
| * @param {number} row - row index of the vector |
| * @return {MatrixRowView} |
| */ |
| rowView(row) { |
| util.checkRowIndex(this, row); |
| return new MLMatrixRowView(this, row); |
| } |
| |
| /** |
| * Returns a view of the column vector with the given index |
| * @param {number} column - column index of the vector |
| * @return {MatrixColumnView} |
| */ |
| columnView(column) { |
| util.checkColumnIndex(this, column); |
| return new MLMatrixColumnView(this, column); |
| } |
| |
| /** |
| * Returns a view of the matrix flipped in the row axis |
| * @return {MatrixFlipRowView} |
| */ |
| flipRowView() { |
| return new MLMatrixFlipRowView(this); |
| } |
| |
| /** |
| * Returns a view of the matrix flipped in the column axis |
| * @return {MatrixFlipColumnView} |
| */ |
| flipColumnView() { |
| return new MLMatrixFlipColumnView(this); |
| } |
| |
| /** |
| * Returns a view of a submatrix giving the index boundaries |
| * @param {number} startRow - first row index of the submatrix |
| * @param {number} endRow - last row index of the submatrix |
| * @param {number} startColumn - first column index of the submatrix |
| * @param {number} endColumn - last column index of the submatrix |
| * @return {MatrixSubView} |
| */ |
| subMatrixView(startRow, endRow, startColumn, endColumn) { |
| return new MLMatrixSubView(this, startRow, endRow, startColumn, endColumn); |
| } |
| |
| /** |
| * Returns a view of the cross of the row indices and the column indices |
| * @example |
| * // resulting vector is [[2], [2]] |
| * var matrix = new Matrix([[1,2,3], [4,5,6]]).selectionView([0, 0], [1]) |
| * @param {Array<number>} rowIndices |
| * @param {Array<number>} columnIndices |
| * @return {MatrixSelectionView} |
| */ |
| selectionView(rowIndices, columnIndices) { |
| return new MLMatrixSelectionView(this, rowIndices, columnIndices); |
| } |
| |
| |
| /** |
| * Calculates and returns the determinant of a matrix as a Number |
| * @example |
| * new Matrix([[1,2,3], [4,5,6]]).det() |
| * @return {number} |
| */ |
| det() { |
| if (this.isSquare()) { |
| var a, b, c, d; |
| if (this.columns === 2) { |
| // 2 x 2 matrix |
| a = this.get(0, 0); |
| b = this.get(0, 1); |
| c = this.get(1, 0); |
| d = this.get(1, 1); |
| |
| return a * d - (b * c); |
| } else if (this.columns === 3) { |
| // 3 x 3 matrix |
| var subMatrix0, subMatrix1, subMatrix2; |
| subMatrix0 = this.selectionView([1, 2], [1, 2]); |
| subMatrix1 = this.selectionView([1, 2], [0, 2]); |
| subMatrix2 = this.selectionView([1, 2], [0, 1]); |
| a = this.get(0, 0); |
| b = this.get(0, 1); |
| c = this.get(0, 2); |
| |
| return a * subMatrix0.det() - b * subMatrix1.det() + c * subMatrix2.det(); |
| } else { |
| // general purpose determinant using the LU decomposition |
| return new LuDecomposition(this).determinant; |
| } |
| |
| } else { |
| throw Error('Determinant can only be calculated for a square matrix.'); |
| } |
| } |
| |
| /** |
| * Returns inverse of a matrix if it exists or the pseudoinverse |
| * @param {number} threshold - threshold for taking inverse of singular values (default = 1e-15) |
| * @return {Matrix} the (pseudo)inverted matrix. |
| */ |
| pseudoInverse(threshold) { |
| if (threshold === undefined) threshold = Number.EPSILON; |
| var svdSolution = new SvDecomposition(this, {autoTranspose: true}); |
| |
| var U = svdSolution.leftSingularVectors; |
| var V = svdSolution.rightSingularVectors; |
| var s = svdSolution.diagonal; |
| |
| for (var i = 0; i < s.length; i++) { |
| if (Math.abs(s[i]) > threshold) { |
| s[i] = 1.0 / s[i]; |
| } else { |
| s[i] = 0.0; |
| } |
| } |
| |
| // convert list to diagonal |
| s = this.constructor[Symbol.species].diag(s); |
| return V.mmul(s.mmul(U.transposeView())); |
| } |
| } |
| |
| Matrix.prototype.klass = 'Matrix'; |
| |
| /** |
| * @private |
| * Check that two matrices have the same dimensions |
| * @param {Matrix} matrix |
| * @param {Matrix} otherMatrix |
| */ |
| function checkDimensions(matrix, otherMatrix) { // eslint-disable-line no-unused-vars |
| if (matrix.rows !== otherMatrix.rows || |
| matrix.columns !== otherMatrix.columns) { |
| throw new RangeError('Matrices dimensions must be equal'); |
| } |
| } |
| |
| function compareNumbers(a, b) { |
| return a - b; |
| } |
| |
| /* |
| Synonyms |
| */ |
| |
| Matrix.random = Matrix.rand; |
| Matrix.diagonal = Matrix.diag; |
| Matrix.prototype.diagonal = Matrix.prototype.diag; |
| Matrix.identity = Matrix.eye; |
| Matrix.prototype.negate = Matrix.prototype.neg; |
| Matrix.prototype.tensorProduct = Matrix.prototype.kroneckerProduct; |
| Matrix.prototype.determinant = Matrix.prototype.det; |
| |
| /* |
| Add dynamically instance and static methods for mathematical operations |
| */ |
| |
| var inplaceOperator = ` |
| (function %name%(value) { |
| if (typeof value === 'number') return this.%name%S(value); |
| return this.%name%M(value); |
| }) |
| `; |
| |
| var inplaceOperatorScalar = ` |
| (function %name%S(value) { |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) %op% value); |
| } |
| } |
| return this; |
| }) |
| `; |
| |
| var inplaceOperatorMatrix = ` |
| (function %name%M(matrix) { |
| matrix = this.constructor.checkMatrix(matrix); |
| checkDimensions(this, matrix); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, this.get(i, j) %op% matrix.get(i, j)); |
| } |
| } |
| return this; |
| }) |
| `; |
| |
| var staticOperator = ` |
| (function %name%(matrix, value) { |
| var newMatrix = new this[Symbol.species](matrix); |
| return newMatrix.%name%(value); |
| }) |
| `; |
| |
| var inplaceMethod = ` |
| (function %name%() { |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, %method%(this.get(i, j))); |
| } |
| } |
| return this; |
| }) |
| `; |
| |
| var staticMethod = ` |
| (function %name%(matrix) { |
| var newMatrix = new this[Symbol.species](matrix); |
| return newMatrix.%name%(); |
| }) |
| `; |
| |
| var inplaceMethodWithArgs = ` |
| (function %name%(%args%) { |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, %method%(this.get(i, j), %args%)); |
| } |
| } |
| return this; |
| }) |
| `; |
| |
| var staticMethodWithArgs = ` |
| (function %name%(matrix, %args%) { |
| var newMatrix = new this[Symbol.species](matrix); |
| return newMatrix.%name%(%args%); |
| }) |
| `; |
| |
| |
| var inplaceMethodWithOneArgScalar = ` |
| (function %name%S(value) { |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, %method%(this.get(i, j), value)); |
| } |
| } |
| return this; |
| }) |
| `; |
| var inplaceMethodWithOneArgMatrix = ` |
| (function %name%M(matrix) { |
| matrix = this.constructor.checkMatrix(matrix); |
| checkDimensions(this, matrix); |
| for (var i = 0; i < this.rows; i++) { |
| for (var j = 0; j < this.columns; j++) { |
| this.set(i, j, %method%(this.get(i, j), matrix.get(i, j))); |
| } |
| } |
| return this; |
| }) |
| `; |
| |
| var inplaceMethodWithOneArg = ` |
| (function %name%(value) { |
| if (typeof value === 'number') return this.%name%S(value); |
| return this.%name%M(value); |
| }) |
| `; |
| |
| var staticMethodWithOneArg = staticMethodWithArgs; |
| |
| var operators = [ |
| // Arithmetic operators |
| ['+', 'add'], |
| ['-', 'sub', 'subtract'], |
| ['*', 'mul', 'multiply'], |
| ['/', 'div', 'divide'], |
| ['%', 'mod', 'modulus'], |
| // Bitwise operators |
| ['&', 'and'], |
| ['|', 'or'], |
| ['^', 'xor'], |
| ['<<', 'leftShift'], |
| ['>>', 'signPropagatingRightShift'], |
| ['>>>', 'rightShift', 'zeroFillRightShift'] |
| ]; |
| |
| var i; |
| |
| for (var operator of operators) { |
| var inplaceOp = eval(fillTemplateFunction(inplaceOperator, {name: operator[1], op: operator[0]})); |
| var inplaceOpS = eval(fillTemplateFunction(inplaceOperatorScalar, {name: operator[1] + 'S', op: operator[0]})); |
| var inplaceOpM = eval(fillTemplateFunction(inplaceOperatorMatrix, {name: operator[1] + 'M', op: operator[0]})); |
| var staticOp = eval(fillTemplateFunction(staticOperator, {name: operator[1]})); |
| for (i = 1; i < operator.length; i++) { |
| Matrix.prototype[operator[i]] = inplaceOp; |
| Matrix.prototype[operator[i] + 'S'] = inplaceOpS; |
| Matrix.prototype[operator[i] + 'M'] = inplaceOpM; |
| Matrix[operator[i]] = staticOp; |
| } |
| } |
| |
| var methods = [ |
| ['~', 'not'] |
| ]; |
| |
| [ |
| 'abs', 'acos', 'acosh', 'asin', 'asinh', 'atan', 'atanh', 'cbrt', 'ceil', |
| 'clz32', 'cos', 'cosh', 'exp', 'expm1', 'floor', 'fround', 'log', 'log1p', |
| 'log10', 'log2', 'round', 'sign', 'sin', 'sinh', 'sqrt', 'tan', 'tanh', 'trunc' |
| ].forEach(function (mathMethod) { |
| methods.push(['Math.' + mathMethod, mathMethod]); |
| }); |
| |
| for (var method of methods) { |
| var inplaceMeth = eval(fillTemplateFunction(inplaceMethod, {name: method[1], method: method[0]})); |
| var staticMeth = eval(fillTemplateFunction(staticMethod, {name: method[1]})); |
| for (i = 1; i < method.length; i++) { |
| Matrix.prototype[method[i]] = inplaceMeth; |
| Matrix[method[i]] = staticMeth; |
| } |
| } |
| |
| var methodsWithArgs = [ |
| ['Math.pow', 1, 'pow'] |
| ]; |
| |
| for (var methodWithArg of methodsWithArgs) { |
| var args = 'arg0'; |
| for (i = 1; i < methodWithArg[1]; i++) { |
| args += `, arg${i}`; |
| } |
| if (methodWithArg[1] !== 1) { |
| var inplaceMethWithArgs = eval(fillTemplateFunction(inplaceMethodWithArgs, { |
| name: methodWithArg[2], |
| method: methodWithArg[0], |
| args: args |
| })); |
| var staticMethWithArgs = eval(fillTemplateFunction(staticMethodWithArgs, {name: methodWithArg[2], args: args})); |
| for (i = 2; i < methodWithArg.length; i++) { |
| Matrix.prototype[methodWithArg[i]] = inplaceMethWithArgs; |
| Matrix[methodWithArg[i]] = staticMethWithArgs; |
| } |
| } else { |
| var tmplVar = { |
| name: methodWithArg[2], |
| args: args, |
| method: methodWithArg[0] |
| }; |
| var inplaceMethod2 = eval(fillTemplateFunction(inplaceMethodWithOneArg, tmplVar)); |
| var inplaceMethodS = eval(fillTemplateFunction(inplaceMethodWithOneArgScalar, tmplVar)); |
| var inplaceMethodM = eval(fillTemplateFunction(inplaceMethodWithOneArgMatrix, tmplVar)); |
| var staticMethod2 = eval(fillTemplateFunction(staticMethodWithOneArg, tmplVar)); |
| for (i = 2; i < methodWithArg.length; i++) { |
| Matrix.prototype[methodWithArg[i]] = inplaceMethod2; |
| Matrix.prototype[methodWithArg[i] + 'M'] = inplaceMethodM; |
| Matrix.prototype[methodWithArg[i] + 'S'] = inplaceMethodS; |
| Matrix[methodWithArg[i]] = staticMethod2; |
| } |
| } |
| } |
| |
| function fillTemplateFunction(template, values) { |
| for (var value in values) { |
| template = template.replace(new RegExp('%' + value + '%', 'g'), values[value]); |
| } |
| return template; |
| } |
| |
| return Matrix; |
| } |
| } |
| |
| |
| // ml-matrix src/views/base |
| let MLMatrixBaseView; |
| { |
| let abstractMatrix = MLMatrixAbstractMatrix; |
| let Matrix = MLMatrixMatrix; |
| |
| class BaseView extends abstractMatrix() { |
| constructor(matrix, rows, columns) { |
| super(); |
| this.matrix = matrix; |
| this.rows = rows; |
| this.columns = columns; |
| } |
| |
| static get [Symbol.species]() { |
| return Matrix.Matrix; |
| } |
| } |
| |
| MLMatrixBaseView = BaseView; |
| } |
| |
| |
| // ml-matrix src/views/column.js |
| let MLMatrixColumnView; |
| { |
| let BaseView = MLMatrixBaseView; |
| |
| class MatrixColumnView extends BaseView { |
| constructor(matrix, column) { |
| super(matrix, matrix.rows, 1); |
| this.column = column; |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(rowIndex, this.column, value); |
| return this; |
| } |
| |
| get(rowIndex) { |
| return this.matrix.get(rowIndex, this.column); |
| } |
| } |
| |
| MLMatrixColumnView = MatrixColumnView; |
| } |
| |
| |
| // ml-matrix src/views/flipColumn.js |
| let MLMatrixFlipColumnView; |
| { |
| let BaseView = MLMatrixBaseView |
| |
| class MatrixFlipColumnView extends BaseView { |
| constructor(matrix) { |
| super(matrix, matrix.rows, matrix.columns); |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(rowIndex, this.columns - columnIndex - 1, value); |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this.matrix.get(rowIndex, this.columns - columnIndex - 1); |
| } |
| } |
| |
| MLMatrixFlipColumnView = MatrixFlipColumnView; |
| } |
| |
| |
| // ml-matrix src/views/flipRow.js |
| let MLMatrixFlipRowView; |
| { |
| let BaseView = MLMatrixBaseView |
| |
| class MatrixFlipRowView extends BaseView { |
| constructor(matrix) { |
| super(matrix, matrix.rows, matrix.columns); |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(this.rows - rowIndex - 1, columnIndex, value); |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this.matrix.get(this.rows - rowIndex - 1, columnIndex); |
| } |
| } |
| |
| MLMatrixFlipRowView = MatrixFlipRowView; |
| } |
| |
| // ml-matrix src/views/row.js |
| let MLMatrixRowView; |
| { |
| let BaseView = MLMatrixBaseView; |
| |
| class MatrixRowView extends BaseView { |
| constructor(matrix, row) { |
| super(matrix, 1, matrix.columns); |
| this.row = row; |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(this.row, columnIndex, value); |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this.matrix.get(this.row, columnIndex); |
| } |
| } |
| |
| MLMatrixRowView = MatrixRowView; |
| } |
| |
| |
| // ml-matrix src/views/selection.js |
| let MLMatrixSelectionView; |
| { |
| let BaseView = MLMatrixBaseView; |
| let util = MLMatrixUtil; |
| |
| class MatrixSelectionView extends BaseView { |
| constructor(matrix, rowIndices, columnIndices) { |
| var indices = util.checkIndices(matrix, rowIndices, columnIndices); |
| super(matrix, indices.row.length, indices.column.length); |
| this.rowIndices = indices.row; |
| this.columnIndices = indices.column; |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(this.rowIndices[rowIndex], this.columnIndices[columnIndex], value); |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this.matrix.get(this.rowIndices[rowIndex], this.columnIndices[columnIndex]); |
| } |
| } |
| |
| MLMatrixSelectionView = MatrixSelectionView; |
| } |
| |
| // ml-matrix src/views/sub.js |
| let MLMatrixSubView; |
| { |
| let BaseView = MLMatrixBaseView; |
| let util = MLMatrixUtil; |
| |
| class MatrixSubView extends BaseView { |
| constructor(matrix, startRow, endRow, startColumn, endColumn) { |
| util.checkRange(matrix, startRow, endRow, startColumn, endColumn); |
| super(matrix, endRow - startRow + 1, endColumn - startColumn + 1); |
| this.startRow = startRow; |
| this.startColumn = startColumn; |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(this.startRow + rowIndex, this.startColumn + columnIndex, value); |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this.matrix.get(this.startRow + rowIndex, this.startColumn + columnIndex); |
| } |
| } |
| |
| MLMatrixSubView = MatrixSubView; |
| } |
| |
| // ml-matrix src/views/transpose.js |
| let MLMatrixTransposeView; |
| { |
| let BaseView = MLMatrixBaseView; |
| |
| class MatrixTransposeView extends BaseView { |
| constructor(matrix) { |
| super(matrix, matrix.columns, matrix.rows); |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this.matrix.set(columnIndex, rowIndex, value); |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this.matrix.get(columnIndex, rowIndex); |
| } |
| } |
| |
| MLMatrixTransposeView = MatrixTransposeView; |
| } |
| |
| // mlmatrix src/matrix.js |
| { |
| let abstractMatrix = MLMatrixAbstractMatrix; |
| let util = MLMatrixUtil; |
| |
| class Matrix extends abstractMatrix(Array) { |
| constructor(nRows, nColumns) { |
| var i; |
| if (arguments.length === 1 && typeof nRows === 'number') { |
| return new Array(nRows); |
| } |
| if (Matrix.isMatrix(nRows)) { |
| return nRows.clone(); |
| } else if (Number.isInteger(nRows) && nRows > 0) { // Create an empty matrix |
| super(nRows); |
| if (Number.isInteger(nColumns) && nColumns > 0) { |
| for (i = 0; i < nRows; i++) { |
| this[i] = new Array(nColumns); |
| } |
| } else { |
| throw new TypeError('nColumns must be a positive integer'); |
| } |
| } else if (Array.isArray(nRows)) { // Copy the values from the 2D array |
| const matrix = nRows; |
| nRows = matrix.length; |
| nColumns = matrix[0].length; |
| if (typeof nColumns !== 'number' || nColumns === 0) { |
| throw new TypeError('Data must be a 2D array with at least one element'); |
| } |
| super(nRows); |
| for (i = 0; i < nRows; i++) { |
| if (matrix[i].length !== nColumns) { |
| throw new RangeError('Inconsistent array dimensions'); |
| } |
| this[i] = [].concat(matrix[i]); |
| } |
| } else { |
| throw new TypeError('First argument must be a positive number or an array'); |
| } |
| this.rows = nRows; |
| this.columns = nColumns; |
| return this; |
| } |
| |
| set(rowIndex, columnIndex, value) { |
| this[rowIndex][columnIndex] = value; |
| return this; |
| } |
| |
| get(rowIndex, columnIndex) { |
| return this[rowIndex][columnIndex]; |
| } |
| |
| /** |
| * Creates an exact and independent copy of the matrix |
| * @return {Matrix} |
| */ |
| clone() { |
| var newMatrix = new this.constructor[Symbol.species](this.rows, this.columns); |
| for (var row = 0; row < this.rows; row++) { |
| for (var column = 0; column < this.columns; column++) { |
| newMatrix.set(row, column, this.get(row, column)); |
| } |
| } |
| return newMatrix; |
| } |
| |
| /** |
| * Removes a row from the given index |
| * @param {number} index - Row index |
| * @return {Matrix} this |
| */ |
| removeRow(index) { |
| util.checkRowIndex(this, index); |
| if (this.rows === 1) { |
| throw new RangeError('A matrix cannot have less than one row'); |
| } |
| this.splice(index, 1); |
| this.rows -= 1; |
| return this; |
| } |
| |
| /** |
| * Adds a row at the given index |
| * @param {number} [index = this.rows] - Row index |
| * @param {Array|Matrix} array - Array or vector |
| * @return {Matrix} this |
| */ |
| addRow(index, array) { |
| if (array === undefined) { |
| array = index; |
| index = this.rows; |
| } |
| util.checkRowIndex(this, index, true); |
| array = util.checkRowVector(this, array, true); |
| this.splice(index, 0, array); |
| this.rows += 1; |
| return this; |
| } |
| |
| /** |
| * Removes a column from the given index |
| * @param {number} index - Column index |
| * @return {Matrix} this |
| */ |
| removeColumn(index) { |
| util.checkColumnIndex(this, index); |
| if (this.columns === 1) { |
| throw new RangeError('A matrix cannot have less than one column'); |
| } |
| for (var i = 0; i < this.rows; i++) { |
| this[i].splice(index, 1); |
| } |
| this.columns -= 1; |
| return this; |
| } |
| |
| /** |
| * Adds a column at the given index |
| * @param {number} [index = this.columns] - Column index |
| * @param {Array|Matrix} array - Array or vector |
| * @return {Matrix} this |
| */ |
| addColumn(index, array) { |
| if (typeof array === 'undefined') { |
| array = index; |
| index = this.columns; |
| } |
| util.checkColumnIndex(this, index, true); |
| array = util.checkColumnVector(this, array); |
| for (var i = 0; i < this.rows; i++) { |
| this[i].splice(index, 0, array[i]); |
| } |
| this.columns += 1; |
| return this; |
| } |
| } |
| |
| MLMatrixMatrix.Matrix = Matrix; |
| Matrix.abstractMatrix = abstractMatrix; |
| } |
| |
| |
| // ml-matrix src/dc/cholesky.js |
| let MLMatrixDCCholesky = {}; |
| { |
| let Matrix = MLMatrixMatrix.Matrix; |
| |
| // https://github.com/lutzroeder/Mapack/blob/master/Source/CholeskyDecomposition.cs |
| function CholeskyDecomposition(value) { |
| if (!(this instanceof CholeskyDecomposition)) { |
| return new CholeskyDecomposition(value); |
| } |
| value = Matrix.checkMatrix(value); |
| if (!value.isSymmetric()) { |
| throw new Error('Matrix is not symmetric'); |
| } |
| |
| var a = value, |
| dimension = a.rows, |
| l = new Matrix(dimension, dimension), |
| positiveDefinite = true, |
| i, j, k; |
| |
| for (j = 0; j < dimension; j++) { |
| var Lrowj = l[j]; |
| var d = 0; |
| for (k = 0; k < j; k++) { |
| var Lrowk = l[k]; |
| var s = 0; |
| for (i = 0; i < k; i++) { |
| s += Lrowk[i] * Lrowj[i]; |
| } |
| Lrowj[k] = s = (a[j][k] - s) / l[k][k]; |
| d = d + s * s; |
| } |
| |
| d = a[j][j] - d; |
| |
| positiveDefinite &= (d > 0); |
| l[j][j] = Math.sqrt(Math.max(d, 0)); |
| for (k = j + 1; k < dimension; k++) { |
| l[j][k] = 0; |
| } |
| } |
| |
| if (!positiveDefinite) { |
| throw new Error('Matrix is not positive definite'); |
| } |
| |
| this.L = l; |
| } |
| |
| CholeskyDecomposition.prototype = { |
| get lowerTriangularMatrix() { |
| return this.L; |
| }, |
| solve: function (value) { |
| value = Matrix.checkMatrix(value); |
| |
| var l = this.L, |
| dimension = l.rows; |
| |
| if (value.rows !== dimension) { |
| throw new Error('Matrix dimensions do not match'); |
| } |
| |
| var count = value.columns, |
| B = value.clone(), |
| i, j, k; |
| |
| for (k = 0; k < dimension; k++) { |
| for (j = 0; j < count; j++) { |
| for (i = 0; i < k; i++) { |
| B[k][j] -= B[i][j] * l[k][i]; |
| } |
| B[k][j] /= l[k][k]; |
| } |
| } |
| |
| for (k = dimension - 1; k >= 0; k--) { |
| for (j = 0; j < count; j++) { |
| for (i = k + 1; i < dimension; i++) { |
| B[k][j] -= B[i][j] * l[i][k]; |
| } |
| B[k][j] /= l[k][k]; |
| } |
| } |
| |
| return B; |
| } |
| }; |
| |
| MLMatrixDCCholesky = CholeskyDecomposition; |
| } |
| |
| |
| // ml-matrix src/dc/evd.js |
| let MLMatrixDCEVD; |
| { |
| const Matrix = MLMatrixMatrix.Matrix; |
| const util = MLMatrixDCUtil; |
| const hypotenuse = util.hypotenuse; |
| const getFilled2DArray = util.getFilled2DArray; |
| |
| const defaultOptions = { |
| assumeSymmetric: false |
| }; |
| |
| // https://github.com/lutzroeder/Mapack/blob/master/Source/EigenvalueDecomposition.cs |
| function EigenvalueDecomposition(matrix, options) { |
| options = Object.assign({}, defaultOptions, options); |
| if (!(this instanceof EigenvalueDecomposition)) { |
| return new EigenvalueDecomposition(matrix, options); |
| } |
| matrix = Matrix.checkMatrix(matrix); |
| if (!matrix.isSquare()) { |
| throw new Error('Matrix is not a square matrix'); |
| } |
| |
| var n = matrix.columns, |
| V = getFilled2DArray(n, n, 0), |
| d = new Array(n), |
| e = new Array(n), |
| value = matrix, |
| i, j; |
| |
| var isSymmetric = false; |
| if (options.assumeSymmetric) { |
| isSymmetric = true; |
| } else { |
| isSymmetric = matrix.isSymmetric(); |
| } |
| |
| if (isSymmetric) { |
| for (i = 0; i < n; i++) { |
| for (j = 0; j < n; j++) { |
| V[i][j] = value.get(i, j); |
| } |
| } |
| tred2(n, e, d, V); |
| tql2(n, e, d, V); |
| } else { |
| var H = getFilled2DArray(n, n, 0), |
| ort = new Array(n); |
| for (j = 0; j < n; j++) { |
| for (i = 0; i < n; i++) { |
| H[i][j] = value.get(i, j); |
| } |
| } |
| orthes(n, H, ort, V); |
| hqr2(n, e, d, V, H); |
| } |
| |
| this.n = n; |
| this.e = e; |
| this.d = d; |
| this.V = V; |
| } |
| |
| EigenvalueDecomposition.prototype = { |
| get realEigenvalues() { |
| return this.d; |
| }, |
| get imaginaryEigenvalues() { |
| return this.e; |
| }, |
| get eigenvectorMatrix() { |
| if (!Matrix.isMatrix(this.V)) { |
| this.V = new Matrix(this.V); |
| } |
| return this.V; |
| }, |
| get diagonalMatrix() { |
| var n = this.n, |
| e = this.e, |
| d = this.d, |
| X = new Matrix(n, n), |
| i, j; |
| for (i = 0; i < n; i++) { |
| for (j = 0; j < n; j++) { |
| X[i][j] = 0; |
| } |
| X[i][i] = d[i]; |
| if (e[i] > 0) { |
| X[i][i + 1] = e[i]; |
| } else if (e[i] < 0) { |
| X[i][i - 1] = e[i]; |
| } |
| } |
| return X; |
| } |
| }; |
| |
| function tred2(n, e, d, V) { |
| |
| var f, g, h, i, j, k, |
| hh, scale; |
| |
| for (j = 0; j < n; j++) { |
| d[j] = V[n - 1][j]; |
| } |
| |
| for (i = n - 1; i > 0; i--) { |
| scale = 0; |
| h = 0; |
| for (k = 0; k < i; k++) { |
| scale = scale + Math.abs(d[k]); |
| } |
| |
| if (scale === 0) { |
| e[i] = d[i - 1]; |
| for (j = 0; j < i; j++) { |
| d[j] = V[i - 1][j]; |
| V[i][j] = 0; |
| V[j][i] = 0; |
| } |
| } else { |
| for (k = 0; k < i; k++) { |
| d[k] /= scale; |
| h += d[k] * d[k]; |
| } |
| |
| f = d[i - 1]; |
| g = Math.sqrt(h); |
| if (f > 0) { |
| g = -g; |
| } |
| |
| e[i] = scale * g; |
| h = h - f * g; |
| d[i - 1] = f - g; |
| for (j = 0; j < i; j++) { |
| e[j] = 0; |
| } |
| |
| for (j = 0; j < i; j++) { |
| f = d[j]; |
| V[j][i] = f; |
| g = e[j] + V[j][j] * f; |
| for (k = j + 1; k <= i - 1; k++) { |
| g += V[k][j] * d[k]; |
| e[k] += V[k][j] * f; |
| } |
| e[j] = g; |
| } |
| |
| f = 0; |
| for (j = 0; j < i; j++) { |
| e[j] /= h; |
| f += e[j] * d[j]; |
| } |
| |
| hh = f / (h + h); |
| for (j = 0; j < i; j++) { |
| e[j] -= hh * d[j]; |
| } |
| |
| for (j = 0; j < i; j++) { |
| f = d[j]; |
| g = e[j]; |
| for (k = j; k <= i - 1; k++) { |
| V[k][j] -= (f * e[k] + g * d[k]); |
| } |
| d[j] = V[i - 1][j]; |
| V[i][j] = 0; |
| } |
| } |
| d[i] = h; |
| } |
| |
| for (i = 0; i < n - 1; i++) { |
| V[n - 1][i] = V[i][i]; |
| V[i][i] = 1; |
| h = d[i + 1]; |
| if (h !== 0) { |
| for (k = 0; k <= i; k++) { |
| d[k] = V[k][i + 1] / h; |
| } |
| |
| for (j = 0; j <= i; j++) { |
| g = 0; |
| for (k = 0; k <= i; k++) { |
| g += V[k][i + 1] * V[k][j]; |
| } |
| for (k = 0; k <= i; k++) { |
| V[k][j] -= g * d[k]; |
| } |
| } |
| } |
| |
| for (k = 0; k <= i; k++) { |
| V[k][i + 1] = 0; |
| } |
| } |
| |
| for (j = 0; j < n; j++) { |
| d[j] = V[n - 1][j]; |
| V[n - 1][j] = 0; |
| } |
| |
| V[n - 1][n - 1] = 1; |
| e[0] = 0; |
| } |
| |
| function tql2(n, e, d, V) { |
| |
| var g, h, i, j, k, l, m, p, r, |
| dl1, c, c2, c3, el1, s, s2, |
| iter; |
| |
| for (i = 1; i < n; i++) { |
| e[i - 1] = e[i]; |
| } |
| |
| e[n - 1] = 0; |
| |
| var f = 0, |
| tst1 = 0, |
| eps = Math.pow(2, -52); |
| |
| for (l = 0; l < n; l++) { |
| tst1 = Math.max(tst1, Math.abs(d[l]) + Math.abs(e[l])); |
| m = l; |
| while (m < n) { |
| if (Math.abs(e[m]) <= eps * tst1) { |
| break; |
| } |
| m++; |
| } |
| |
| if (m > l) { |
| iter = 0; |
| do { |
| iter = iter + 1; |
| |
| g = d[l]; |
| p = (d[l + 1] - g) / (2 * e[l]); |
| r = hypotenuse(p, 1); |
| if (p < 0) { |
| r = -r; |
| } |
| |
| d[l] = e[l] / (p + r); |
| d[l + 1] = e[l] * (p + r); |
| dl1 = d[l + 1]; |
| h = g - d[l]; |
| for (i = l + 2; i < n; i++) { |
| d[i] -= h; |
| } |
| |
| f = f + h; |
| |
| p = d[m]; |
| c = 1; |
| c2 = c; |
| c3 = c; |
| el1 = e[l + 1]; |
| s = 0; |
| s2 = 0; |
| for (i = m - 1; i >= l; i--) { |
| c3 = c2; |
| c2 = c; |
| s2 = s; |
| g = c * e[i]; |
| h = c * p; |
| r = hypotenuse(p, e[i]); |
| e[i + 1] = s * r; |
| s = e[i] / r; |
| c = p / r; |
| p = c * d[i] - s * g; |
| d[i + 1] = h + s * (c * g + s * d[i]); |
| |
| for (k = 0; k < n; k++) { |
| h = V[k][i + 1]; |
| V[k][i + 1] = s * V[k][i] + c * h; |
| V[k][i] = c * V[k][i] - s * h; |
| } |
| } |
| |
| p = -s * s2 * c3 * el1 * e[l] / dl1; |
| e[l] = s * p; |
| d[l] = c * p; |
| |
| } |
| while (Math.abs(e[l]) > eps * tst1); |
| } |
| d[l] = d[l] + f; |
| e[l] = 0; |
| } |
| |
| for (i = 0; i < n - 1; i++) { |
| k = i; |
| p = d[i]; |
| for (j = i + 1; j < n; j++) { |
| if (d[j] < p) { |
| k = j; |
| p = d[j]; |
| } |
| } |
| |
| if (k !== i) { |
| d[k] = d[i]; |
| d[i] = p; |
| for (j = 0; j < n; j++) { |
| p = V[j][i]; |
| V[j][i] = V[j][k]; |
| V[j][k] = p; |
| } |
| } |
| } |
| } |
| |
| function orthes(n, H, ort, V) { |
| |
| var low = 0, |
| high = n - 1, |
| f, g, h, i, j, m, |
| scale; |
| |
| for (m = low + 1; m <= high - 1; m++) { |
| scale = 0; |
| for (i = m; i <= high; i++) { |
| scale = scale + Math.abs(H[i][m - 1]); |
| } |
| |
| if (scale !== 0) { |
| h = 0; |
| for (i = high; i >= m; i--) { |
| ort[i] = H[i][m - 1] / scale; |
| h += ort[i] * ort[i]; |
| } |
| |
| g = Math.sqrt(h); |
| if (ort[m] > 0) { |
| g = -g; |
| } |
| |
| h = h - ort[m] * g; |
| ort[m] = ort[m] - g; |
| |
| for (j = m; j < n; j++) { |
| f = 0; |
| for (i = high; i >= m; i--) { |
| f += ort[i] * H[i][j]; |
| } |
| |
| f = f / h; |
| for (i = m; i <= high; i++) { |
| H[i][j] -= f * ort[i]; |
| } |
| } |
| |
| for (i = 0; i <= high; i++) { |
| f = 0; |
| for (j = high; j >= m; j--) { |
| f += ort[j] * H[i][j]; |
| } |
| |
| f = f / h; |
| for (j = m; j <= high; j++) { |
| H[i][j] -= f * ort[j]; |
| } |
| } |
| |
| ort[m] = scale * ort[m]; |
| H[m][m - 1] = scale * g; |
| } |
| } |
| |
| for (i = 0; i < n; i++) { |
| for (j = 0; j < n; j++) { |
| V[i][j] = (i === j ? 1 : 0); |
| } |
| } |
| |
| for (m = high - 1; m >= low + 1; m--) { |
| if (H[m][m - 1] !== 0) { |
| for (i = m + 1; i <= high; i++) { |
| ort[i] = H[i][m - 1]; |
| } |
| |
| for (j = m; j <= high; j++) { |
| g = 0; |
| for (i = m; i <= high; i++) { |
| g += ort[i] * V[i][j]; |
| } |
| |
| g = (g / ort[m]) / H[m][m - 1]; |
| for (i = m; i <= high; i++) { |
| V[i][j] += g * ort[i]; |
| } |
| } |
| } |
| } |
| } |
| |
| function hqr2(nn, e, d, V, H) { |
| var n = nn - 1, |
| low = 0, |
| high = nn - 1, |
| eps = Math.pow(2, -52), |
| exshift = 0, |
| norm = 0, |
| p = 0, |
| q = 0, |
| r = 0, |
| s = 0, |
| z = 0, |
| iter = 0, |
| i, j, k, l, m, t, w, x, y, |
| ra, sa, vr, vi, |
| notlast, cdivres; |
| |
| for (i = 0; i < nn; i++) { |
| if (i < low || i > high) { |
| d[i] = H[i][i]; |
| e[i] = 0; |
| } |
| |
| for (j = Math.max(i - 1, 0); j < nn; j++) { |
| norm = norm + Math.abs(H[i][j]); |
| } |
| } |
| |
| while (n >= low) { |
| l = n; |
| while (l > low) { |
| s = Math.abs(H[l - 1][l - 1]) + Math.abs(H[l][l]); |
| if (s === 0) { |
| s = norm; |
| } |
| if (Math.abs(H[l][l - 1]) < eps * s) { |
| break; |
| } |
| l--; |
| } |
| |
| if (l === n) { |
| H[n][n] = H[n][n] + exshift; |
| d[n] = H[n][n]; |
| e[n] = 0; |
| n--; |
| iter = 0; |
| } else if (l === n - 1) { |
| w = H[n][n - 1] * H[n - 1][n]; |
| p = (H[n - 1][n - 1] - H[n][n]) / 2; |
| q = p * p + w; |
| z = Math.sqrt(Math.abs(q)); |
| H[n][n] = H[n][n] + exshift; |
| H[n - 1][n - 1] = H[n - 1][n - 1] + exshift; |
| x = H[n][n]; |
| |
| if (q >= 0) { |
| z = (p >= 0) ? (p + z) : (p - z); |
| d[n - 1] = x + z; |
| d[n] = d[n - 1]; |
| if (z !== 0) { |
| d[n] = x - w / z; |
| } |
| e[n - 1] = 0; |
| e[n] = 0; |
| x = H[n][n - 1]; |
| s = Math.abs(x) + Math.abs(z); |
| p = x / s; |
| q = z / s; |
| r = Math.sqrt(p * p + q * q); |
| p = p / r; |
| q = q / r; |
| |
| for (j = n - 1; j < nn; j++) { |
| z = H[n - 1][j]; |
| H[n - 1][j] = q * z + p * H[n][j]; |
| H[n][j] = q * H[n][j] - p * z; |
| } |
| |
| for (i = 0; i <= n; i++) { |
| z = H[i][n - 1]; |
| H[i][n - 1] = q * z + p * H[i][n]; |
| H[i][n] = q * H[i][n] - p * z; |
| } |
| |
| for (i = low; i <= high; i++) { |
| z = V[i][n - 1]; |
| V[i][n - 1] = q * z + p * V[i][n]; |
| V[i][n] = q * V[i][n] - p * z; |
| } |
| } else { |
| d[n - 1] = x + p; |
| d[n] = x + p; |
| e[n - 1] = z; |
| e[n] = -z; |
| } |
| |
| n = n - 2; |
| iter = 0; |
| } else { |
| x = H[n][n]; |
| y = 0; |
| w = 0; |
| if (l < n) { |
| y = H[n - 1][n - 1]; |
| w = H[n][n - 1] * H[n - 1][n]; |
| } |
| |
| if (iter === 10) { |
| exshift += x; |
| for (i = low; i <= n; i++) { |
| H[i][i] -= x; |
| } |
| s = Math.abs(H[n][n - 1]) + Math.abs(H[n - 1][n - 2]); |
| x = y = 0.75 * s; |
| w = -0.4375 * s * s; |
| } |
| |
| if (iter === 30) { |
| s = (y - x) / 2; |
| s = s * s + w; |
| if (s > 0) { |
| s = Math.sqrt(s); |
| if (y < x) { |
| s = -s; |
| } |
| s = x - w / ((y - x) / 2 + s); |
| for (i = low; i <= n; i++) { |
| H[i][i] -= s; |
| } |
| exshift += s; |
| x = y = w = 0.964; |
| } |
| } |
| |
| iter = iter + 1; |
| |
| m = n - 2; |
| while (m >= l) { |
| z = H[m][m]; |
| r = x - z; |
| s = y - z; |
| p = (r * s - w) / H[m + 1][m] + H[m][m + 1]; |
| q = H[m + 1][m + 1] - z - r - s; |
| r = H[m + 2][m + 1]; |
| s = Math.abs(p) + Math.abs(q) + Math.abs(r); |
| p = p / s; |
| q = q / s; |
| r = r / s; |
| if (m === l) { |
| break; |
| } |
| if (Math.abs(H[m][m - 1]) * (Math.abs(q) + Math.abs(r)) < eps * (Math.abs(p) * (Math.abs(H[m - 1][m - 1]) + Math.abs(z) + Math.abs(H[m + 1][m + 1])))) { |
| break; |
| } |
| m--; |
| } |
| |
| for (i = m + 2; i <= n; i++) { |
| H[i][i - 2] = 0; |
| if (i > m + 2) { |
| H[i][i - 3] = 0; |
| } |
| } |
| |
| for (k = m; k <= n - 1; k++) { |
| notlast = (k !== n - 1); |
| if (k !== m) { |
| p = H[k][k - 1]; |
| q = H[k + 1][k - 1]; |
| r = (notlast ? H[k + 2][k - 1] : 0); |
| x = Math.abs(p) + Math.abs(q) + Math.abs(r); |
| if (x !== 0) { |
| p = p / x; |
| q = q / x; |
| r = r / x; |
| } |
| } |
| |
| if (x === 0) { |
| break; |
| } |
| |
| s = Math.sqrt(p * p + q * q + r * r); |
| if (p < 0) { |
| s = -s; |
| } |
| |
| if (s !== 0) { |
| if (k !== m) { |
| H[k][k - 1] = -s * x; |
| } else if (l !== m) { |
| H[k][k - 1] = -H[k][k - 1]; |
| } |
| |
| p = p + s; |
| x = p / s; |
| y = q / s; |
| z = r / s; |
| q = q / p; |
| r = r / p; |
| |
| for (j = k; j < nn; j++) { |
| p = H[k][j] + q * H[k + 1][j]; |
| if (notlast) { |
| p = p + r * H[k + 2][j]; |
| H[k + 2][j] = H[k + 2][j] - p * z; |
| } |
| |
| H[k][j] = H[k][j] - p * x; |
| H[k + 1][j] = H[k + 1][j] - p * y; |
| } |
| |
| for (i = 0; i <= Math.min(n, k + 3); i++) { |
| p = x * H[i][k] + y * H[i][k + 1]; |
| if (notlast) { |
| p = p + z * H[i][k + 2]; |
| H[i][k + 2] = H[i][k + 2] - p * r; |
| } |
| |
| H[i][k] = H[i][k] - p; |
| H[i][k + 1] = H[i][k + 1] - p * q; |
| } |
| |
| for (i = low; i <= high; i++) { |
| p = x * V[i][k] + y * V[i][k + 1]; |
| if (notlast) { |
| p = p + z * V[i][k + 2]; |
| V[i][k + 2] = V[i][k + 2] - p * r; |
| } |
| |
| V[i][k] = V[i][k] - p; |
| V[i][k + 1] = V[i][k + 1] - p * q; |
| } |
| } |
| } |
| } |
| } |
| |
| if (norm === 0) { |
| return; |
| } |
| |
| for (n = nn - 1; n >= 0; n--) { |
| p = d[n]; |
| q = e[n]; |
| |
| if (q === 0) { |
| l = n; |
| H[n][n] = 1; |
| for (i = n - 1; i >= 0; i--) { |
| w = H[i][i] - p; |
| r = 0; |
| for (j = l; j <= n; j++) { |
| r = r + H[i][j] * H[j][n]; |
| } |
| |
| if (e[i] < 0) { |
| z = w; |
| s = r; |
| } else { |
| l = i; |
| if (e[i] === 0) { |
| H[i][n] = (w !== 0) ? (-r / w) : (-r / (eps * norm)); |
| } else { |
| x = H[i][i + 1]; |
| y = H[i + 1][i]; |
| q = (d[i] - p) * (d[i] - p) + e[i] * e[i]; |
| t = (x * s - z * r) / q; |
| H[i][n] = t; |
| H[i + 1][n] = (Math.abs(x) > Math.abs(z)) ? ((-r - w * t) / x) : ((-s - y * t) / z); |
| } |
| |
| t = Math.abs(H[i][n]); |
| if ((eps * t) * t > 1) { |
| for (j = i; j <= n; j++) { |
| H[j][n] = H[j][n] / t; |
| } |
| } |
| } |
| } |
| } else if (q < 0) { |
| l = n - 1; |
| |
| if (Math.abs(H[n][n - 1]) > Math.abs(H[n - 1][n])) { |
| H[n - 1][n - 1] = q / H[n][n - 1]; |
| H[n - 1][n] = -(H[n][n] - p) / H[n][n - 1]; |
| } else { |
| cdivres = cdiv(0, -H[n - 1][n], H[n - 1][n - 1] - p, q); |
| H[n - 1][n - 1] = cdivres[0]; |
| H[n - 1][n] = cdivres[1]; |
| } |
| |
| H[n][n - 1] = 0; |
| H[n][n] = 1; |
| for (i = n - 2; i >= 0; i--) { |
| ra = 0; |
| sa = 0; |
| for (j = l; j <= n; j++) { |
| ra = ra + H[i][j] * H[j][n - 1]; |
| sa = sa + H[i][j] * H[j][n]; |
| } |
| |
| w = H[i][i] - p; |
| |
| if (e[i] < 0) { |
| z = w; |
| r = ra; |
| s = sa; |
| } else { |
| l = i; |
| if (e[i] === 0) { |
| cdivres = cdiv(-ra, -sa, w, q); |
| H[i][n - 1] = cdivres[0]; |
| H[i][n] = cdivres[1]; |
| } else { |
| x = H[i][i + 1]; |
| y = H[i + 1][i]; |
| vr = (d[i] - p) * (d[i] - p) + e[i] * e[i] - q * q; |
| vi = (d[i] - p) * 2 * q; |
| if (vr === 0 && vi === 0) { |
| vr = eps * norm * (Math.abs(w) + Math.abs(q) + Math.abs(x) + Math.abs(y) + Math.abs(z)); |
| } |
| cdivres = cdiv(x * r - z * ra + q * sa, x * s - z * sa - q * ra, vr, vi); |
| H[i][n - 1] = cdivres[0]; |
| H[i][n] = cdivres[1]; |
| if (Math.abs(x) > (Math.abs(z) + Math.abs(q))) { |
| H[i + 1][n - 1] = (-ra - w * H[i][n - 1] + q * H[i][n]) / x; |
| H[i + 1][n] = (-sa - w * H[i][n] - q * H[i][n - 1]) / x; |
| } else { |
| cdivres = cdiv(-r - y * H[i][n - 1], -s - y * H[i][n], z, q); |
| H[i + 1][n - 1] = cdivres[0]; |
| H[i + 1][n] = cdivres[1]; |
| } |
| } |
| |
| t = Math.max(Math.abs(H[i][n - 1]), Math.abs(H[i][n])); |
| if ((eps * t) * t > 1) { |
| for (j = i; j <= n; j++) { |
| H[j][n - 1] = H[j][n - 1] / t; |
| H[j][n] = H[j][n] / t; |
| } |
| } |
| } |
| } |
| } |
| } |
| |
| for (i = 0; i < nn; i++) { |
| if (i < low || i > high) { |
| for (j = i; j < nn; j++) { |
| V[i][j] = H[i][j]; |
| } |
| } |
| } |
| |
| for (j = nn - 1; j >= low; j--) { |
| for (i = low; i <= high; i++) { |
| z = 0; |
| for (k = low; k <= Math.min(j, high); k++) { |
| z = z + V[i][k] * H[k][j]; |
| } |
| V[i][j] = z; |
| } |
| } |
| } |
| |
| function cdiv(xr, xi, yr, yi) { |
| var r, d; |
| if (Math.abs(yr) > Math.abs(yi)) { |
| r = yi / yr; |
| d = yr + r * yi; |
| return [(xr + r * xi) / d, (xi - r * xr) / d]; |
| } else { |
| r = yr / yi; |
| d = yi + r * yr; |
| return [(r * xr + xi) / d, (r * xi - xr) / d]; |
| } |
| } |
| |
| MLMatrixDCEVD = EigenvalueDecomposition; |
| } |
| |
| |
| // ml-matrix src/dc/qr.js |
| let MLMatrixDCQR; |
| { |
| let Matrix = MLMatrixMatrix.Matrix; |
| let hypotenuse = MLMatrixDCUtil.hypotenuse; |
| |
| //https://github.com/lutzroeder/Mapack/blob/master/Source/QrDecomposition.cs |
| function QrDecomposition(value) { |
| if (!(this instanceof QrDecomposition)) { |
| return new QrDecomposition(value); |
| } |
| value = Matrix.checkMatrix(value); |
| |
| var qr = value.clone(), |
| m = value.rows, |
| n = value.columns, |
| rdiag = new Array(n), |
| i, j, k, s; |
| |
| for (k = 0; k < n; k++) { |
| var nrm = 0; |
| for (i = k; i < m; i++) { |
| nrm = hypotenuse(nrm, qr[i][k]); |
| } |
| if (nrm !== 0) { |
| if (qr[k][k] < 0) { |
| nrm = -nrm; |
| } |
| for (i = k; i < m; i++) { |
| qr[i][k] /= nrm; |
| } |
| qr[k][k] += 1; |
| for (j = k + 1; j < n; j++) { |
| s = 0; |
| for (i = k; i < m; i++) { |
| s += qr[i][k] * qr[i][j]; |
| } |
| s = -s / qr[k][k]; |
| for (i = k; i < m; i++) { |
| qr[i][j] += s * qr[i][k]; |
| } |
| } |
| } |
| rdiag[k] = -nrm; |
| } |
| |
| this.QR = qr; |
| this.Rdiag = rdiag; |
| } |
| |
| QrDecomposition.prototype = { |
| solve: function (value) { |
| value = Matrix.checkMatrix(value); |
| |
| var qr = this.QR, |
| m = qr.rows; |
| |
| if (value.rows !== m) { |
| throw new Error('Matrix row dimensions must agree'); |
| } |
| if (!this.isFullRank()) { |
| throw new Error('Matrix is rank deficient'); |
| } |
| |
| var count = value.columns; |
| var X = value.clone(); |
| var n = qr.columns; |
| var i, j, k, s; |
| |
| for (k = 0; k < n; k++) { |
| for (j = 0; j < count; j++) { |
| s = 0; |
| for (i = k; i < m; i++) { |
| s += qr[i][k] * X[i][j]; |
| } |
| s = -s / qr[k][k]; |
| for (i = k; i < m; i++) { |
| X[i][j] += s * qr[i][k]; |
| } |
| } |
| } |
| for (k = n - 1; k >= 0; k--) { |
| for (j = 0; j < count; j++) { |
| X[k][j] /= this.Rdiag[k]; |
| } |
| for (i = 0; i < k; i++) { |
| for (j = 0; j < count; j++) { |
| X[i][j] -= X[k][j] * qr[i][k]; |
| } |
| } |
| } |
| |
| return X.subMatrix(0, n - 1, 0, count - 1); |
| }, |
| isFullRank: function () { |
| var columns = this.QR.columns; |
| for (var i = 0; i < columns; i++) { |
| if (this.Rdiag[i] === 0) { |
| return false; |
| } |
| } |
| return true; |
| }, |
| get upperTriangularMatrix() { |
| var qr = this.QR, |
| n = qr.columns, |
| X = new Matrix(n, n), |
| i, j; |
| for (i = 0; i < n; i++) { |
| for (j = 0; j < n; j++) { |
| if (i < j) { |
| X[i][j] = qr[i][j]; |
| } else if (i === j) { |
| X[i][j] = this.Rdiag[i]; |
| } else { |
| X[i][j] = 0; |
| } |
| } |
| } |
| return X; |
| }, |
| get orthogonalMatrix() { |
| var qr = this.QR, |
| rows = qr.rows, |
| columns = qr.columns, |
| X = new Matrix(rows, columns), |
| i, j, k, s; |
| |
| for (k = columns - 1; k >= 0; k--) { |
| for (i = 0; i < rows; i++) { |
| X[i][k] = 0; |
| } |
| X[k][k] = 1; |
| for (j = k; j < columns; j++) { |
| if (qr[k][k] !== 0) { |
| s = 0; |
| for (i = k; i < rows; i++) { |
| s += qr[i][k] * X[i][j]; |
| } |
| |
| s = -s / qr[k][k]; |
| |
| for (i = k; i < rows; i++) { |
| X[i][j] += s * qr[i][k]; |
| } |
| } |
| } |
| } |
| return X; |
| } |
| }; |
| |
| MLMatrixDCQR = QrDecomposition; |
| } |
| |
| // ml-matric src/decompositions.js |
| let MLMatrixDecompositions = {}; |
| { |
| let Matrix = MLMatrixMatrix.Matrix; |
| |
| let SingularValueDecomposition = MLMatrixDCSVD; |
| let EigenvalueDecomposition = MLMatrixDCEVD; |
| let LuDecomposition = MLMatrixDCLU; |
| let QrDecomposition = MLMatrixDCQR |
| let CholeskyDecomposition = MLMatrixDCCholesky; |
| |
| function inverse(matrix) { |
| matrix = Matrix.checkMatrix(matrix); |
| return solve(matrix, Matrix.eye(matrix.rows)); |
| } |
| |
| /** |
| * Returns the inverse |
| * @memberOf Matrix |
| * @static |
| * @param {Matrix} matrix |
| * @return {Matrix} matrix |
| * @alias inv |
| */ |
| Matrix.inverse = Matrix.inv = inverse; |
| |
| /** |
| * Returns the inverse |
| * @memberOf Matrix |
| * @static |
| * @param {Matrix} matrix |
| * @return {Matrix} matrix |
| * @alias inv |
| */ |
| Matrix.prototype.inverse = Matrix.prototype.inv = function () { |
| return inverse(this); |
| }; |
| |
| function solve(leftHandSide, rightHandSide) { |
| leftHandSide = Matrix.checkMatrix(leftHandSide); |
| rightHandSide = Matrix.checkMatrix(rightHandSide); |
| return leftHandSide.isSquare() ? new LuDecomposition(leftHandSide).solve(rightHandSide) : new QrDecomposition(leftHandSide).solve(rightHandSide); |
| } |
| |
| Matrix.solve = solve; |
| Matrix.prototype.solve = function (other) { |
| return solve(this, other); |
| }; |
| |
| MLMatrixDecompositions = { |
| SingularValueDecomposition: SingularValueDecomposition, |
| SVD: SingularValueDecomposition, |
| EigenvalueDecomposition: EigenvalueDecomposition, |
| EVD: EigenvalueDecomposition, |
| LuDecomposition: LuDecomposition, |
| LU: LuDecomposition, |
| QrDecomposition: QrDecomposition, |
| QR: QrDecomposition, |
| CholeskyDecomposition: CholeskyDecomposition, |
| CHO: CholeskyDecomposition, |
| inverse: inverse, |
| solve: solve |
| }; |
| } |
| |
| // ml-matrix src/index.js |
| let MLMatrix = {}; |
| { |
| MLMatrix = MLMatrixMatrix.Matrix; |
| MLMatrix.Decompositions = MLMatrix.DC = MLMatrixDecompositions; |
| } |
| |
| // feedforward-neural-networks utils.js |
| let FeedforwardNeuralNetworksUtils; |
| { |
| let Matrix = MLMatrix; |
| |
| /** |
| * @private |
| * Retrieves the sum at each row of the given matrix. |
| * @param {Matrix} matrix |
| * @return {Matrix} |
| */ |
| function sumRow(matrix) { |
| var sum = Matrix.zeros(matrix.rows, 1); |
| for (var i = 0; i < matrix.rows; ++i) { |
| for (var j = 0; j < matrix.columns; ++j) { |
| sum[i][0] += matrix[i][j]; |
| } |
| } |
| return sum; |
| } |
| |
| /** |
| * @private |
| * Retrieves the sum at each column of the given matrix. |
| * @param {Matrix} matrix |
| * @return {Matrix} |
| */ |
| function sumCol(matrix) { |
| var sum = Matrix.zeros(1, matrix.columns); |
| for (var i = 0; i < matrix.rows; ++i) { |
| for (var j = 0; j < matrix.columns; ++j) { |
| sum[0][j] += matrix[i][j]; |
| } |
| } |
| return sum; |
| } |
| |
| /** |
| * @private |
| * Method that given an array of labels(predictions), returns two dictionaries, one to transform from labels to |
| * numbers and other in the reverse way |
| * @param {Array} array |
| * @return {object} |
| */ |
| function dictOutputs(array) { |
| var inputs = {}, outputs = {}, l = array.length, index = 0; |
| for (var i = 0; i < l; i += 1) { |
| if (inputs[array[i]] === undefined) { |
| inputs[array[i]] = index; |
| outputs[index] = array[i]; |
| index++; |
| } |
| } |
| |
| return { |
| inputs: inputs, |
| outputs: outputs |
| }; |
| } |
| |
| FeedforwardNeuralNetworksUtils = { |
| dictOutputs: dictOutputs, |
| sumCol: sumCol, |
| sumRow: sumRow |
| }; |
| } |
| |
| // feedforward-neural-networks activationFunctions.js |
| let FeedforwardNeuralNetworksActivationFunctions; |
| { |
| function logistic(val) { |
| return 1 / (1 + Math.exp(-val)); |
| } |
| |
| function expELU(val, param) { |
| return val < 0 ? param * (Math.exp(val) - 1) : val; |
| } |
| |
| function softExponential(val, param) { |
| if (param < 0) { |
| return -Math.log(1 - param * (val + param)) / param; |
| } |
| if (param > 0) { |
| return ((Math.exp(param * val) - 1) / param) + param; |
| } |
| return val; |
| } |
| |
| function softExponentialPrime(val, param) { |
| if (param < 0) { |
| return 1 / (1 - param * (param + val)); |
| } else { |
| return Math.exp(param * val); |
| } |
| } |
| |
| const ACTIVATION_FUNCTIONS = { |
| 'tanh': { |
| activation: Math.tanh, |
| derivate: val => 1 - (val * val) |
| }, |
| 'identity': { |
| activation: val => val, |
| derivate: () => 1 |
| }, |
| 'logistic': { |
| activation: logistic, |
| derivate: val => logistic(val) * (1 - logistic(val)) |
| }, |
| 'arctan': { |
| activation: Math.atan, |
| derivate: val => 1 / (val * val + 1) |
| }, |
| 'softsign': { |
| activation: val => val / (1 + Math.abs(val)), |
| derivate: val => 1 / ((1 + Math.abs(val)) * (1 + Math.abs(val))) |
| }, |
| 'relu': { |
| activation: val => val < 0 ? 0 : val, |
| derivate: val => val < 0 ? 0 : 1 |
| }, |
| 'softplus': { |
| activation: val => Math.log(1 + Math.exp(val)), |
| derivate: val => 1 / (1 + Math.exp(-val)) |
| }, |
| 'bent': { |
| activation: val => ((Math.sqrt(val * val + 1) - 1) / 2) + val, |
| derivate: val => (val / (2 * Math.sqrt(val * val + 1))) + 1 |
| }, |
| 'sinusoid': { |
| activation: Math.sin, |
| derivate: Math.cos |
| }, |
| 'sinc': { |
| activation: val => val === 0 ? 1 : Math.sin(val) / val, |
| derivate: val => val === 0 ? 0 : (Math.cos(val) / val) - (Math.sin(val) / (val * val)) |
| }, |
| 'gaussian': { |
| activation: val => Math.exp(-(val * val)), |
| derivate: val => -2 * val * Math.exp(-(val * val)) |
| }, |
| 'parametric-relu': { |
| activation: (val, param) => val < 0 ? param * val : val, |
| derivate: (val, param) => val < 0 ? param : 1 |
| }, |
| 'exponential-elu': { |
| activation: expELU, |
| derivate: (val, param) => val < 0 ? expELU(val, param) + param : 1 |
| }, |
| 'soft-exponential': { |
| activation: softExponential, |
| derivate: softExponentialPrime |
| } |
| }; |
| |
| FeedforwardNeuralNetworksActivationFunctions = ACTIVATION_FUNCTIONS; |
| } |
| |
| // feedforward-neural-networks Layer.js |
| let FeedforwardNeuralNetworksLayer; |
| { |
| let Matrix = MLMatrix; |
| |
| let Utils = FeedforwardNeuralNetworksUtils; |
| const ACTIVATION_FUNCTIONS = FeedforwardNeuralNetworksActivationFunctions; |
| |
| class Layer { |
| /** |
| * @private |
| * Create a new layer with the given options |
| * @param {object} options |
| * @param {number} [options.inputSize] - Number of conections that enter the neurons. |
| * @param {number} [options.outputSize] - Number of conections that leave the neurons. |
| * @param {number} [options.regularization] - Regularization parameter. |
| * @param {number} [options.epsilon] - Learning rate parameter. |
| * @param {string} [options.activation] - Activation function parameter from the FeedForwardNeuralNetwork class. |
| * @param {number} [options.activationParam] - Activation parameter if needed. |
| */ |
| constructor(options) { |
| this.inputSize = options.inputSize; |
| this.outputSize = options.outputSize; |
| this.regularization = options.regularization; |
| this.epsilon = options.epsilon; |
| this.activation = options.activation; |
| this.activationParam = options.activationParam; |
| |
| var selectedFunction = ACTIVATION_FUNCTIONS[options.activation]; |
| var params = selectedFunction.activation.length; |
| |
| var actFunction = params > 1 ? val => selectedFunction.activation(val, options.activationParam) : selectedFunction.activation; |
| var derFunction = params > 1 ? val => selectedFunction.derivate(val, options.activationParam) : selectedFunction.derivate; |
| |
| this.activationFunction = function (i, j) { |
| this[i][j] = actFunction(this[i][j]); |
| }; |
| this.derivate = function (i, j) { |
| this[i][j] = derFunction(this[i][j]); |
| }; |
| |
| if (options.model) { |
| // load model |
| this.W = Matrix.checkMatrix(options.W); |
| this.b = Matrix.checkMatrix(options.b); |
| |
| } else { |
| // default constructor |
| |
| this.W = Matrix.rand(this.inputSize, this.outputSize); |
| this.b = Matrix.zeros(1, this.outputSize); |
| |
| this.W.apply(function (i, j) { |
| this[i][j] /= Math.sqrt(options.inputSize); |
| }); |
| } |
| } |
| |
| /** |
| * @private |
| * propagate the given input through the current layer. |
| * @param {Matrix} X - input. |
| * @return {Matrix} output at the current layer. |
| */ |
| forward(X) { |
| var z = X.mmul(this.W).addRowVector(this.b); |
| z.apply(this.activationFunction); |
| this.a = z.clone(); |
| return z; |
| } |
| |
| /** |
| * @private |
| * apply backpropagation algorithm at the current layer |
| * @param {Matrix} delta - delta values estimated at the following layer. |
| * @param {Matrix} a - 'a' values from the following layer. |
| * @return {Matrix} the new delta values for the next layer. |
| */ |
| backpropagation(delta, a) { |
| this.dW = a.transposeView().mmul(delta); |
| this.db = Utils.sumCol(delta); |
| |
| var aCopy = a.clone(); |
| return delta.mmul(this.W.transposeView()).mul(aCopy.apply(this.derivate)); |
| } |
| |
| /** |
| * @private |
| * Function that updates the weights at the current layer with the derivatives. |
| */ |
| update() { |
| this.dW.add(this.W.clone().mul(this.regularization)); |
| this.W.add(this.dW.mul(-this.epsilon)); |
| this.b.add(this.db.mul(-this.epsilon)); |
| } |
| |
| /** |
| * @private |
| * Export the current layer to JSON. |
| * @return {object} model |
| */ |
| toJSON() { |
| return { |
| model: 'Layer', |
| inputSize: this.inputSize, |
| outputSize: this.outputSize, |
| regularization: this.regularization, |
| epsilon: this.epsilon, |
| activation: this.activation, |
| W: this.W, |
| b: this.b |
| }; |
| } |
| |
| /** |
| * @private |
| * Creates a new Layer with the given model. |
| * @param {object} model |
| * @return {Layer} |
| */ |
| static load(model) { |
| if (model.model !== 'Layer') { |
| throw new RangeError('the current model is not a Layer model'); |
| } |
| return new Layer(model); |
| } |
| |
| } |
| |
| FeedforwardNeuralNetworksLayer = Layer; |
| } |
| |
| // feedforward-neural-networks OutputLayer.js |
| let FeedforwardNeuralNetworksOutputLayer; |
| { |
| let Layer = FeedforwardNeuralNetworksLayer; |
| |
| class OutputLayer extends Layer { |
| constructor(options) { |
| super(options); |
| |
| this.activationFunction = function (i, j) { |
| this[i][j] = Math.exp(this[i][j]); |
| }; |
| } |
| |
| static load(model) { |
| if (model.model !== 'Layer') { |
| throw new RangeError('the current model is not a Layer model'); |
| } |
| |
| return new OutputLayer(model); |
| } |
| } |
| |
| FeedforwardNeuralNetworksOutputLayer = OutputLayer; |
| } |
| |
| // feedforward-neural-networks FeedForwardNeuralNetwork.js |
| let FeedforwardNeuralNetwork; |
| { |
| const Matrix = MLMatrix; |
| |
| const Layer = FeedforwardNeuralNetworksLayer; |
| const OutputLayer = FeedforwardNeuralNetworksOutputLayer; |
| const Utils = FeedforwardNeuralNetworksUtils; |
| const ACTIVATION_FUNCTIONS = FeedforwardNeuralNetworksActivationFunctions; |
| |
| class FeedForwardNeuralNetworks { |
| |
| /** |
| * Create a new Feedforword neural network model. |
| * @param {object} options |
| * @param {Array} [options.hiddenLayers=[10]] - Array that contains the sizes of the hidden layers. |
| * @oaram {number} [options.iterations=50] - Number of iterations at the training step. |
| * @param {number} [options.learningRate=0.01] - Learning rate of the neural net (also known as epsilon). |
| * @poram {number} [options.regularization=0.01] - Regularization parameter af the neural net. |
| * @poram {string} [options.activation='tanh'] - activation function to be used. (options: 'tanh'(default), |
| * 'identity', 'logistic', 'arctan', 'softsign', 'relu', 'softplus', 'bent', 'sinusoid', 'sinc', 'gaussian'). |
| * (single-parametric options: 'parametric-relu', 'exponential-relu', 'soft-exponential'). |
| * @param {number} [options.activationParam=1] - if the selected activation function needs a parameter. |
| */ |
| constructor(options) { |
| options = options || {}; |
| if (options.model) { |
| // load network |
| this.hiddenLayers = options.hiddenLayers; |
| this.iterations = options.iterations; |
| this.learningRate = options.learningRate; |
| this.regularization = options.regularization; |
| this.dicts = options.dicts; |
| this.activation = options.activation; |
| this.activationParam = options.activationParam; |
| this.model = new Array(options.layers.length); |
| |
| for (var i = 0; i < this.model.length - 1; ++i) { |
| this.model[i] = Layer.load(options.layers[i]); |
| } |
| this.model[this.model.length - 1] = OutputLayer.load(options.layers[this.model.length - 1]); |
| } else { |
| // default constructor |
| this.hiddenLayers = options.hiddenLayers === undefined ? [10] : options.hiddenLayers; |
| this.iterations = options.iterations === undefined ? 50 : options.iterations; |
| |
| this.learningRate = options.learningRate === undefined ? 0.01 : options.learningRate; |
| //this.momentum = options.momentum === undefined ? 0.1 : options.momentum; |
| this.regularization = options.regularization === undefined ? 0.01 : options.regularization; |
| |
| this.activation = options.activation === undefined ? 'tanh' : options.activation; |
| this.activationParam = options.activationParam === undefined ? 1 : options.activationParam; |
| if (!(this.activation in Object.keys(ACTIVATION_FUNCTIONS))) { |
| this.activation = 'tanh'; |
| } |
| } |
| } |
| |
| /** |
| * @private |
| * Function that build and initialize the neural net. |
| * @param {number} inputSize - total of features to fit. |
| * @param {number} outputSize - total of labels of the prediction set. |
| */ |
| buildNetwork(inputSize, outputSize) { |
| var size = 2 + (this.hiddenLayers.length - 1); |
| this.model = new Array(size); |
| |
| // input layer |
| this.model[0] = new Layer({ |
| inputSize: inputSize, |
| outputSize: this.hiddenLayers[0], |
| activation: this.activation, |
| activationParam: this.activationParam, |
| regularization: this.regularization, |
| epsilon: this.learningRate |
| }); |
| |
| // hidden layers |
| for (var i = 1; i < this.hiddenLayers.length; ++i) { |
| this.model[i] = new Layer({ |
| inputSize: this.hiddenLayers[i - 1], |
| outputSize: this.hiddenLayers[i], |
| activation: this.activation, |
| activationParam: this.activationParam, |
| regularization: this.regularization, |
| epsilon: this.learningRate |
| }); |
| } |
| |
| // output layer |
| this.model[size - 1] = new OutputLayer({ |
| inputSize: this.hiddenLayers[this.hiddenLayers.length - 1], |
| outputSize: outputSize, |
| activation: this.activation, |
| activationParam: this.activationParam, |
| regularization: this.regularization, |
| epsilon: this.learningRate |
| }); |
| } |
| |
| /** |
| * Train the neural net with the given features and labels. |
| * @param {Matrix|Array} features |
| * @param {Matrix|Array} labels |
| */ |
| train(features, labels) { |
| features = Matrix.checkMatrix(features); |
| this.dicts = Utils.dictOutputs(labels); |
| |
| var inputSize = features.columns; |
| var outputSize = Object.keys(this.dicts.inputs).length; |
| |
| this.buildNetwork(inputSize, outputSize); |
| |
| for (var i = 0; i < this.iterations; ++i) { |
| var probabilities = this.propagate(features); |
| this.backpropagation(features, labels, probabilities); |
| } |
| } |
| |
| /** |
| * @private |
| * Propagate the input(training set) and retrives the probabilities of each class. |
| * @param {Matrix} X |
| * @return {Matrix} probabilities of each class. |
| */ |
| propagate(X) { |
| var input = X; |
| for (var i = 0; i < this.model.length; ++i) { |
| //console.log(i); |
| input = this.model[i].forward(input); |
| } |
| |
| // get probabilities |
| return input.divColumnVector(Utils.sumRow(input)); |
| } |
| |
| /** |
| * @private |
| * Function that applies the backpropagation algorithm on each layer of the network |
| * in order to fit the features and labels. |
| * @param {Matrix} features |
| * @param {Array} labels |
| * @param {Matrix} probabilities - probabilities of each class of the feature set. |
| */ |
| backpropagation(features, labels, probabilities) { |
| for (var i = 0; i < probabilities.length; ++i) { |
| probabilities[i][this.dicts.inputs[labels[i]]] -= 1; |
| } |
| |
| // remember, the last delta doesn't matter |
| var delta = probabilities; |
| for (i = this.model.length - 1; i >= 0; --i) { |
| var a = i > 0 ? this.model[i - 1].a : features; |
| delta = this.model[i].backpropagation(delta, a); |
| } |
| |
| for (i = 0; i < this.model.length; ++i) { |
| this.model[i].update(); |
| } |
| } |
| |
| /** |
| * Predict the output given the feature set. |
| * @param {Array|Matrix} features |
| * @return {Array} |
| */ |
| predict(features) { |
| features = Matrix.checkMatrix(features); |
| var outputs = new Array(features.rows); |
| var probabilities = this.propagate(features); |
| for (var i = 0; i < features.rows; ++i) { |
| outputs[i] = this.dicts.outputs[probabilities.maxRowIndex(i)[1]]; |
| } |
| |
| return outputs; |
| } |
| |
| /** |
| * Export the current model to JSOM. |
| * @return {object} model |
| */ |
| toJSON() { |
| var model = { |
| model: 'FNN', |
| hiddenLayers: this.hiddenLayers, |
| iterations: this.iterations, |
| learningRate: this.learningRate, |
| regularization: this.regularization, |
| activation: this.activation, |
| activationParam: this.activationParam, |
| dicts: this.dicts, |
| layers: new Array(this.model.length) |
| }; |
| |
| for (var i = 0; i < this.model.length; ++i) { |
| model.layers[i] = this.model[i].toJSON(); |
| } |
| |
| return model; |
| } |
| |
| /** |
| * Load a Feedforward Neural Network with the current model. |
| * @param {object} model |
| * @return {FeedForwardNeuralNetworks} |
| */ |
| static load(model) { |
| if (model.model !== 'FNN') { |
| throw new RangeError('the current model is not a feed forward network'); |
| } |
| |
| return new FeedForwardNeuralNetworks(model); |
| } |
| } |
| |
| FeedforwardNeuralNetwork = FeedForwardNeuralNetworks; |
| } |