Nnmf

examples/linearlyDependent.js

@@ -0,0 +1,11 @@
+import Matrix from '../src';
+
+const test = [
+  [1, 0, 0, 0, 0, 0, 0, 0, 0, 0],
+  [0, 1, 0, 0, 0, 0, 0, 0, 0, 0],
+  [1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
+  [1, 1, 1, 0, 0, 0, 0, 0, 0, 0],
+  [1, 1, 1, 1, 0, 0, 0, 0, 0, 0],
+  [1, 1, 1, 1, 1, 0, 0, 0, 0, 0],
+  [5, 4, 3, 2, 1, 0, 0, 0, 0, 0]
+];

examples/nnmf.js

@@ -0,0 +1,21 @@
+import Matrix from '../src';
+
+const test = [
+  [1, 0, 0, 0, 0, 0, 0, 0, 0, 0],
+  [1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
+  [1, 1, 1, 0, 0, 0, 0, 0, 0, 0],
+  [1, 1, 1, 1, 0, 0, 0, 0, 0, 0],
+  [1, 1, 1, 1, 1, 0, 0, 0, 0, 0],
+  [1, 1, 1, 1, 1, 1, 0, 0, 0, 0],
+  [1, 1, 1, 1, 1, 1, 1, 0, 0, 0],
+  [1, 1, 1, 1, 1, 1, 1, 1, 0, 0],
+  [1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
+  [1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
+];
+const targets = [
+  [1, 0, 0, 0, 0, 0, 0, 0, 0, 0],
+  [0, 1, 0, 0, 0, 0, 0, 0, 0, 0],
+  [0, 1, 1, 0, 0, 0, 0, 0, 0, 0],
+  [10, 9, 8, 7, 6, 5, 4, 3, 2, 1],
+  [4, 4, 2, 2, 0, 0, 0, 0, 0, 0]
+];

package.json

@@ -57,6 +57,7 @@
     "eslint-config-cheminfo": "^1.18.0",
     "eslint-plugin-import": "^2.14.0",
     "eslint-plugin-jest": "^21.24.1",
+    "esm": "^3.0.84",
     "jest": "^23.6.0",
     "jest-matcher-deep-close-to": "^1.3.0",
     "mathjs": "^5.1.2",

src/__tests__/decompositions/nnmf.js

@@ -0,0 +1,22 @@
+import { toBeDeepCloseTo } from 'jest-matcher-deep-close-to';
+
+import { Matrix, NNMF } from '../..';
+
+expect.extend({ toBeDeepCloseTo });
+
+describe('Non-negative Matrix Factorization', () => {
+  it('factorization', () => {
+    let A = new Matrix([
+      [1, 0, 0, 0, 1],
+      [0, 0, 0, 0, 0],
+      [0, 0, 0, 0, 0],
+      [0, 0, 0, 0, 0],
+      [1, 0, 0, 0, 1]
+    ]);
+
+    let nA = new NNMF(A, 4);
+    expect(nA.X.mmul(nA.Y)).toBeDeepCloseTo(A);
+    expect(nA.error).toBeDeepCloseTo(Matrix.zeros(5, 5));
+    expect(nA.positivity).toEqual(true);
+  });
+});

src/abstractMatrix.js

@@ -132,16 +132,18 @@ export default function AbstractMatrix(superCtor) {
          * @param {number} rows - Number of rows
          * @param {number} columns - Number of columns
          * @param {number} [maxValue=1000] - Maximum value
+         * @param {number} [minValue=0] - Minimum value
          * @param {function} [rng=Math.random] - Random number generator
          * @return {Matrix} The new matrix
          */
-    static randInt(rows, columns, maxValue, rng) {
+    static randInt(rows, columns, maxValue, minValue, rng) {
       if (maxValue === undefined) maxValue = 1000;
+      if (minValue === undefined) minValue = 0;
       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);
+          var value = Math.floor(rng() * (maxValue - minValue) + minValue);
           matrix.set(i, j, value);
         }
       }

src/dc/nnmf.js

@@ -0,0 +1,137 @@
+import { Matrix, WrapperMatrix2D } from '../index';
+
+/**
+ *
+ * @class NNMF
+ * @link http://papers.nips.cc/paper/1861-algorithms-for-non-negative-matrix-factorization.pdf
+ * @param {Matrix} A
+ * @param {number} r
+ * @param {object} [options={}]
+ * @param {number} [options.numberIterations=1000]
+ */
+export default class NNMF {
+  constructor(A, r, options = {}) {
+    const { numberIterations = 1000 } = options;
+    A = WrapperMatrix2D.checkMatrix(A);
+    var m = A.rows;
+    var n = A.columns;
+
+    this.r = r;
+    this.A = A;
+    this.m = m;
+    this.n = n;
+    this.X = Matrix.rand(m, r);
+    this.Y = Matrix.rand(r, n);
+
+    doNnmf.call(this, numberIterations);
+  }
+
+  /**
+   * Compute the error
+   */
+  get error() {
+    let error = new Matrix(this.m, this.n);
+    let A2 = this.X.mmul(this.Y);
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.n; j++) {
+        error.set(
+          i,
+          j,
+          Math.abs(A2.get(i, j) - this.A.get(i, j)) /
+            (this.A.get(i, j) !== 0 ? this.A.get(i, j) : Number.EPSILON)
+        );
+      }
+    }
+    return (error);
+  }
+  get positivity() {
+    let positive = true;
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.r; j++) {
+        if (this.X.get(i, j) < 0) {
+          positive = false;
+        }
+      }
+    }
+    for (let i = 0; i < this.r; i++) {
+      for (let j = 0; j < this.n; j++) {
+        if (this.X.get(i, j) < 0) {
+          positive = false;
+        }
+      }
+    }
+    return (positive);
+  }
+}
+
+/**
+ * Do the NNMF of a matrix A into two matrix X and Y
+ * @param {number} [numberIterations=1000]
+ */
+function doNnmf(numberIterations = 1000) {
+  let A2 = this.X.mmul(this.Y);
+  let numX = Matrix.empty(this.m, this.r);
+  let denumX = Matrix.empty(this.m, this.r);
+
+  let numY = Matrix.empty(this.r, this.n);
+  let denumY = Matrix.empty(this.r, this.n);
+
+  let temp1 = Matrix.empty(this.m, this.n);
+  let temp2 = Matrix.empty(this.n, this.m);
+
+  for (let a = 0; a < numberIterations; a++) {
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.n; j++) {
+        temp1.set(i, j, Math.pow(A2.get(i, j), -2) * this.A.get(i, j));
+        temp2.set(i, j, Math.pow(A2.get(i, j), -1));
+      }
+    }
+    numX = temp1.mmul(this.Y.transpose());
+    denumX = temp2.mmul(this.Y.transpose());
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.r; j++) {
+        numX.set(i, j, numX.get(i, j) + Number.EPSILON);
+        denumX.set(i, j, denumX.get(i, j) + Number.EPSILON);
+      }
+    }
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.r; j++) {
+        this.X.set(
+          i,
+          j,
+          (this.X.get(i, j) * numX.get(i, j)) / denumX.get(i, j)
+        );
+      }
+    }
+    A2 = this.X.mmul(this.Y);
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.n; j++) {
+        A2.set(i, j, A2.get(i, j) + Number.EPSILON);
+      }
+    }
+
+    for (let i = 0; i < this.m; i++) {
+      for (let j = 0; j < this.n; j++) {
+        temp1.set(i, j, Math.pow(A2.get(i, j), -2) * this.A.get(i, j));
+        temp2.set(i, j, Math.pow(A2.get(i, j), -1));
+      }
+    }
+    numY = this.X.transpose().mmul(temp1);
+    denumY = this.X.transpose().mmul(temp2);
+    for (let i = 0; i < this.r; i++) {
+      for (let j = 0; j < this.n; j++) {
+        numY.set(i, j, numY.get(i, j) + Number.EPSILON);
+        denumY.set(i, j, denumY.get(i, j) + Number.EPSILON);
+      }
+    }
+    for (let i = 0; i < this.r; i++) {
+      for (let j = 0; j < this.n; j++) {
+        this.Y.set(
+          i,
+          j,
+          (this.Y.get(i, j) * numY.get(i, j)) / denumY.get(i, j)
+        );
+      }
+    }
+  }
+}

src/index.js

@@ -21,3 +21,4 @@ export {
 } from './dc/cholesky.js';
 export { default as LuDecomposition, default as LU } from './dc/lu.js';
 export { default as QrDecomposition, default as QR } from './dc/qr.js';
+export { default as NNMFDecomposition, default as NNMF } from './dc/nnmf.js';