OPTICSをJavaScriptで実装した

はじめに

「色々な機械学習処理をブラウザ上で試せるサイトを作った」中で実装したモデルの解説の六回目です。

今回はOPTICSの実装について解説します。

デモはこちらから。（TaskをClusteringにして、ModelのOPTICSを選択）
実際のコードはoptics.jsにあります。

なお、可視化部分については一切触れません。

概説

OPTICSについては、Wikipediaの擬似コードをそのままコードに起こしただけとなります。

ただし、内部で優先度付きキューを使用する必要があるため、その実装も行っています。

優先度付きキュー

優先度付きキューの実装はシンプルで、以下の関数を用意します。

• push
  データを追加する。引数は値と優先度の二つ。
• move
  データの優先度を変更する。
• shift
  優先度の高いデータを取得する。

注意すべき点としては、優先度に距離を直接設定しているため、「優先度の値が小さい方が優先度が高い」ということです。

内部実装は、値と優先度のペアの配列を保持するようにしています。

コード

優先度付きキュー

class PriorityQueue {
    constructor(arr) {
        this._value = arr || []
    }

    get length() {
        return this._value.length;
    }

    _sort() {
        this._value.sort((a, b) => a[1] - b[1])
    }

    push(value, priority) {
        this._value.push([value, priority])
        this._sort()
    }

    move(value, priority) {
        for (let i = 0; i < this.length; i++) {
            if (this._value[i][0] === value) {
                this._value[i][1] = priority;
                this._sort()
                return
            }
        }
        this.push(value, priority)
    }

    shift() {
        const [value, priority] = this._value.shift();
        return value;
    }
}

OPTICS

class OPTICS {
    // https://en.wikipedia.org/wiki/OPTICS_algorithm
    constructor(eps = Infinity, minPts = 5, metric = 'euclid') {
        this._eps = eps;
        this._minPts = minPts;

        this._metric = metric
        switch (this._metric) {
        case 'euclid':
            this._d = (a, b) => Math.sqrt(a.reduce((s, v, i) => s + (v - b[i]) ** 2, 0));
            break
        case 'manhattan':
            this._d = (a, b) => a.reduce((s, v, i) => s + Math.abs(v - b[i]), 0)
            break
        case 'chebyshev':
            this._d = (a, b) => Math.max(...a.map((v, i) => Math.abs(v - b[i])))
            break;
        }
    }

    fit(datas) {
        const n = datas.length;
        const d = Array(n);
        for (let i = 0; i < n; d[i++] = Array(n));
        for (let i = 0; i < n; i++) {
            for (let j = 0; j < i; j++) {
                const v = this._d(datas[i], datas[j]);
                d[i][j] = d[j][i] = v;
            }
        }
        const getNeighbors = (i) => {
            const neighbors = [];
            for (let k = 0; k < n; k++) {
                if (d[i][k] < this._eps) neighbors.push(k);
            }
            return neighbors
        }
        const coreDist = (i) => {
            const neighbors = getNeighbors(i).map(k => d[i][k]);
            if (neighbors.length <= this._minPts) return null;
            neighbors.sort((a, b) => a - b);
            return neighbors[this._minPts];
        }
        const reachabilityDist = (o, i) => {
            const cd = coreDist(i);
            if (cd === null) return cd;
            return Math.max(cd, d[i][o])
        }

        const processed = Array(n).fill(false);
        const rd = Array(n).fill(null);
        const update = (n, p, seeds) => {
            const cd = coreDist(p);
            if (cd === null) return
            for (const o of n) {
                if (processed[o]) continue
                const nrd = Math.max(cd, d[p][o]);
                if (rd[o] === null) {
                    rd[o] = nrd;
                    seeds.push(o, nrd);
                } else if (nrd < rd[o]) {
                    rd[o] = nrd;
                    seeds.move(o, nrd)
                }
            }
        }

        this._core_distance = [];
        for (let p = 0; p < n; p++) {
            if (processed[p]) continue;
            const neighbors = getNeighbors(p);
            processed[p] = true;
            const cd = coreDist(p)
            this._core_distance.push([p, cd])
            if (cd !== null) {
                const seeds = new PriorityQueue()
                update(neighbors, p, seeds)
                while (seeds.length > 0) {
                    const q = seeds.shift();
                    const nd = getNeighbors(q);
                    processed[q] = true
                    const cdq = coreDist(q)
                    this._core_distance.push([q, cdq])
                    if (cdq !== null) {
                        update(nd, q, seeds);
                    }
                }
            }
        }
    }

    predict(threshold = 0.1) {
        let c = 0;
        const n = this._core_distance.length
        const clusters = Array(n)
        for (let i = 0; i < n; i++) {
            const [k, d] = this._core_distance[i]
            clusters[k] = d > threshold ? ++c : c;
        }
        return clusters
    }
}