แก้ coding style
This commit is contained in:
Pitchaya Boonsarngsuk
310
src/barnesHut.js
310
src/barnesHut.js
@@ -1,158 +1,152 @@
|
||||
import constant from "./constant";
|
||||
import jiggle from "./jiggle";
|
||||
import {x, y} from "./xy";
|
||||
import {quadtree} from "d3-quadtree";
|
||||
|
||||
/**
|
||||
* The refinement of the existing Barnes-Hut implementation in D3
|
||||
* to fit the use case of the project. Previously the algorithm stored
|
||||
* strength as internal node, now the random child is stored as internal
|
||||
* node and the force calculations are done between the node and that internal
|
||||
* object if they are sufficiently far away.
|
||||
* The check to see if the nodes are far away was also changed to the one described in original Barnes-Hut paper.
|
||||
* @return {force} calculated forces.
|
||||
*/
|
||||
export default function() {
|
||||
var nodes,
|
||||
node,
|
||||
alpha,
|
||||
distance = constant(300),
|
||||
theta = 0.5;
|
||||
|
||||
/**
|
||||
* Constructs a quadtree at every iteration and apply the forces by visiting
|
||||
* each node in a tree.
|
||||
* @param {number} _ - controls the stopping of the
|
||||
* particle simulations.
|
||||
*/
|
||||
function force(_) {
|
||||
var i, n = nodes.length, tree = quadtree(nodes, x, y).visitAfter(accumulate);
|
||||
for (alpha = _, i = 0; i < n; ++i) {
|
||||
node = nodes[i], tree.visit(apply);
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Function used during the tree construction to fill out the nodes with
|
||||
* correct data. Internal nodes acquire the random child while the leaf
|
||||
* nodes accumulate forces from coincident quadrants.
|
||||
* @param {quadrant} quad - node representing the quadrant in quadtree.
|
||||
*/
|
||||
function accumulate(quad) {
|
||||
var q, d, children = [];
|
||||
|
||||
// For internal nodes, accumulate forces from child quadrants.
|
||||
if (quad.length) {
|
||||
for (var i = 0; i < 4; ++i) {
|
||||
if ((q = quad[i]) && (d = q.data)) {
|
||||
children.push(d);
|
||||
}
|
||||
}
|
||||
// Choose a random child.
|
||||
quad.data = children[Math.floor(Math.random() * children.length)];
|
||||
quad.x = quad.data.x;
|
||||
quad.y = quad.data.y;
|
||||
}
|
||||
|
||||
// For leaf nodes, accumulate forces from coincident quadrants.
|
||||
else {
|
||||
q = quad;
|
||||
q.x = q.data.x;
|
||||
q.y = q.data.y;
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Function that applies the forces for each node. If the objects are
|
||||
* far away, the approximation is made. Otherwise, forces are calculated
|
||||
* directly between the nodes.
|
||||
* @param {quadrant} quad - node representing the quadrant in quadtree.
|
||||
* @param {number} x1 - lower x bound of the node.
|
||||
* @param {number} _ - lower y bound of the node.
|
||||
* @param {number} x2 - upper x bound of the node.
|
||||
* @return {boolean} - true if the approximation was applied.
|
||||
*/
|
||||
function apply(quad, x1, _, x2) {
|
||||
|
||||
var x = quad.data.x + quad.data.vx - node.x - node.vx,
|
||||
y = quad.data.y + quad.data.vy - node.y - node.vy,
|
||||
w = x2 - x1,
|
||||
l = Math.sqrt(x * x + y * y);
|
||||
|
||||
// Apply the Barnes-Hut approximation if possible.
|
||||
// Limit forces for very close nodes; randomize direction if coincident.
|
||||
if (w / l < theta) {
|
||||
if (x === 0) x = jiggle(), l += x * x;
|
||||
if (y === 0) y = jiggle(), l += y * y;
|
||||
if (quad.data) {
|
||||
l = (l - +distance(node, quad.data)) / l * alpha;
|
||||
x *= l, y *= l;
|
||||
quad.data.vx -= x;
|
||||
quad.data.vy -= y;
|
||||
node.vx += x;
|
||||
node.vy += y;
|
||||
}
|
||||
return true;
|
||||
}
|
||||
|
||||
// Otherwise, process points directly.
|
||||
else if (quad.length) return;
|
||||
|
||||
// Limit forces for very close nodes; randomize direction if coincident.
|
||||
if (quad.data !== node || quad.next) {
|
||||
if (x === 0) x = jiggle(), l += x * x;
|
||||
if (y === 0) y = jiggle(), l += y * y;
|
||||
}
|
||||
|
||||
do if (quad.data !== node) {
|
||||
l = (l - +distance(node, quad.data)) / l * alpha;
|
||||
x *= l, y *= l;
|
||||
quad.data.vx -= x;
|
||||
quad.data.vy -= y;
|
||||
node.vx += x;
|
||||
node.vy += y;
|
||||
} while (quad = quad.next);
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculates the stress. Basically, it computes the difference between
|
||||
* high dimensional distance and real distance. The lower the stress is,
|
||||
* the better layout.
|
||||
* @return {number} - stress of the layout.
|
||||
*/
|
||||
function getStress() {
|
||||
var totalDiffSq = 0, totalHighDistSq = 0;
|
||||
for (var i = 0, source, target, realDist, highDist; i < nodes.length; i++) {
|
||||
for (var j = 0; j < nodes.length; j++) {
|
||||
if (i !== j) {
|
||||
source = nodes[i], target = nodes[j];
|
||||
realDist = Math.hypot(target.x-source.x, target.y-source.y);
|
||||
highDist = +distance(nodes[i], nodes[j]);
|
||||
totalDiffSq += Math.pow(realDist-highDist, 2);
|
||||
totalHighDistSq += highDist * highDist;
|
||||
}
|
||||
}
|
||||
}
|
||||
return Math.sqrt(totalDiffSq/totalHighDistSq);
|
||||
}
|
||||
|
||||
// API for initializing the algorithm, setting parameters and querying
|
||||
// metrics.
|
||||
force.initialize = function(_) {
|
||||
nodes = _;
|
||||
};
|
||||
|
||||
force.distance = function(_) {
|
||||
return arguments.length ? (distance = typeof _ === "function" ? _ : constant(+_), force) : distance;
|
||||
};
|
||||
|
||||
force.theta = function(_) {
|
||||
return arguments.length ? (theta = _, force) : theta;
|
||||
};
|
||||
|
||||
force.stress = function() {
|
||||
return getStress();
|
||||
};
|
||||
|
||||
return force;
|
||||
}
|
||||
import constant from "./constant";
|
||||
import jiggle from "./jiggle";
|
||||
import {x, y} from "./xy";
|
||||
import {quadtree} from "d3-quadtree";
|
||||
|
||||
/**
|
||||
* The refinement of the existing Barnes-Hut implementation in D3
|
||||
* to fit the use case of the project. Previously the algorithm stored
|
||||
* strength as internal node, now the random child is stored as internal
|
||||
* node and the force calculations are done between the node and that internal
|
||||
* object if they are sufficiently far away.
|
||||
* The check to see if the nodes are far away was also changed to the one described in original Barnes-Hut paper.
|
||||
* @return {force} calculated forces.
|
||||
*/
|
||||
export default function() {
|
||||
var nodes,
|
||||
node,
|
||||
alpha,
|
||||
distance = constant(300),
|
||||
theta = 0.5;
|
||||
|
||||
/**
|
||||
* Constructs a quadtree at every iteration and apply the forces by visiting
|
||||
* each node in a tree.
|
||||
* @param {number} _ - controls the stopping of the
|
||||
* particle simulations.
|
||||
*/
|
||||
function force(_) {
|
||||
var i, n = nodes.length, tree = quadtree(nodes, x, y).visitAfter(accumulate);
|
||||
for (alpha = _, i = 0; i < n; ++i) {
|
||||
node = nodes[i], tree.visit(apply);
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Function used during the tree construction to fill out the nodes with
|
||||
* correct data. Internal nodes acquire the random child while the leaf
|
||||
* nodes accumulate forces from coincident quadrants.
|
||||
* @param {quadrant} quad - node representing the quadrant in quadtree.
|
||||
*/
|
||||
function accumulate(quad) {
|
||||
var q, d, children = [];
|
||||
|
||||
// For internal nodes, accumulate forces from child quadrants.
|
||||
if (quad.length) {
|
||||
for (var i = 0; i < 4; ++i) {
|
||||
if ((q = quad[i]) && (d = q.data)) {
|
||||
children.push(d);
|
||||
}
|
||||
}
|
||||
// Choose a random child.
|
||||
quad.data = children[Math.floor(Math.random() * children.length)];
|
||||
quad.x = quad.data.x;
|
||||
quad.y = quad.data.y;
|
||||
} else { // For leaf nodes, accumulate forces from coincident quadrants.
|
||||
q = quad;
|
||||
q.x = q.data.x;
|
||||
q.y = q.data.y;
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Function that applies the forces for each node. If the objects are
|
||||
* far away, the approximation is made. Otherwise, forces are calculated
|
||||
* directly between the nodes.
|
||||
* @param {quadrant} quad - node representing the quadrant in quadtree.
|
||||
* @param {number} x1 - lower x bound of the node.
|
||||
* @param {number} _ - lower y bound of the node.
|
||||
* @param {number} x2 - upper x bound of the node.
|
||||
* @return {boolean} - true if the approximation was applied.
|
||||
*/
|
||||
function apply(quad, x1, _, x2) {
|
||||
|
||||
var x = quad.data.x + quad.data.vx - node.x - node.vx,
|
||||
y = quad.data.y + quad.data.vy - node.y - node.vy,
|
||||
w = x2 - x1,
|
||||
l = Math.sqrt(x * x + y * y);
|
||||
|
||||
// Apply the Barnes-Hut approximation if possible.
|
||||
// Limit forces for very close nodes; randomize direction if coincident.
|
||||
if (w / l < theta) {
|
||||
if (x === 0) x = jiggle(), l += x * x;
|
||||
if (y === 0) y = jiggle(), l += y * y;
|
||||
if (quad.data) {
|
||||
l = (l - +distance(node, quad.data)) / l * alpha;
|
||||
x *= l, y *= l;
|
||||
quad.data.vx -= x;
|
||||
quad.data.vy -= y;
|
||||
node.vx += x;
|
||||
node.vy += y;
|
||||
}
|
||||
return true;
|
||||
} else if (quad.length) return; // Otherwise, process points directly.
|
||||
|
||||
// Limit forces for very close nodes; randomize direction if coincident.
|
||||
if (quad.data !== node || quad.next) {
|
||||
if (x === 0) x = jiggle(), l += x * x;
|
||||
if (y === 0) y = jiggle(), l += y * y;
|
||||
}
|
||||
|
||||
do if (quad.data !== node) {
|
||||
l = (l - +distance(node, quad.data)) / l * alpha;
|
||||
x *= l, y *= l;
|
||||
quad.data.vx -= x;
|
||||
quad.data.vy -= y;
|
||||
node.vx += x;
|
||||
node.vy += y;
|
||||
} while (quad = quad.next);
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculates the stress. Basically, it computes the difference between
|
||||
* high dimensional distance and real distance. The lower the stress is,
|
||||
* the better layout.
|
||||
* @return {number} - stress of the layout.
|
||||
*/
|
||||
function getStress() {
|
||||
var totalDiffSq = 0, totalHighDistSq = 0;
|
||||
for (var i = 0, source, target, realDist, highDist; i < nodes.length; i++) {
|
||||
for (var j = 0; j < nodes.length; j++) {
|
||||
if (i !== j) {
|
||||
source = nodes[i], target = nodes[j];
|
||||
realDist = Math.hypot(target.x-source.x, target.y-source.y);
|
||||
highDist = +distance(nodes[i], nodes[j]);
|
||||
totalDiffSq += Math.pow(realDist-highDist, 2);
|
||||
totalHighDistSq += highDist * highDist;
|
||||
}
|
||||
}
|
||||
}
|
||||
return Math.sqrt(totalDiffSq/totalHighDistSq);
|
||||
}
|
||||
|
||||
// API for initializing the algorithm, setting parameters and querying
|
||||
// metrics.
|
||||
force.initialize = function(_) {
|
||||
nodes = _;
|
||||
};
|
||||
|
||||
force.distance = function(_) {
|
||||
return arguments.length ? (distance = typeof _ === "function" ? _ : constant(+_), force) : distance;
|
||||
};
|
||||
|
||||
force.theta = function(_) {
|
||||
return arguments.length ? (theta = _, force) : theta;
|
||||
};
|
||||
|
||||
force.stress = function() {
|
||||
return getStress();
|
||||
};
|
||||
|
||||
return force;
|
||||
}
|
||||
|
||||
@@ -1,203 +1,203 @@
|
||||
import {dispatch} from "d3-dispatch";
|
||||
import constant from "./constant";
|
||||
import interpBruteForce from "./interpolation/interpBruteForce";
|
||||
import interpolationPivots from "./interpolation/interpolationPivots";
|
||||
import {takeSampleFrom} from "./interpolation/helpers";
|
||||
|
||||
/**
|
||||
* An implementation of Chalmers, Morrison, and Ross' 2002 hybrid layout
|
||||
* algorithm with an option to use the 2003 pivot-based near neighbour searching
|
||||
* method.
|
||||
* It performs the 1996 neighbour sampling spring simulation model, on a
|
||||
* "sample set", sqrt(n) samples of the data.
|
||||
* Other data points are then interpolated into the model.
|
||||
* Finally, another spring simulation may be performed on the entire dataset to
|
||||
* clean up the model.
|
||||
* @param {object} sim - D3 Simulation object
|
||||
* @param {object} forceS - Pre-configured D3 force object for the sample set.
|
||||
The ending handler will be re-configured.
|
||||
Neighbour sampling force is expected, but other D3
|
||||
forces may also work.
|
||||
* @param {object} forceF - Pre-configured D3 force object for the simultion ran
|
||||
on the entire dataset at the end.
|
||||
Neighbour sampling force is expected, but other D3
|
||||
forces may also work.
|
||||
The force should not have any ending condition.
|
||||
*/
|
||||
export default function (sim, forceS, forceF) {
|
||||
var
|
||||
SAMPLE_ITERATIONS = 300,
|
||||
FULL_ITERATIONS = 20,
|
||||
interpDistanceFn,
|
||||
NUM_PIVOTS = 0,
|
||||
INTERP_FINE_ITS = 20,
|
||||
sample = [],
|
||||
remainder = [],
|
||||
simulation = sim,
|
||||
forceSample = forceS,
|
||||
forceFull = forceF,
|
||||
event = d3.dispatch("sampleTick", "fullTick", "startInterp", "end"),
|
||||
initAlready = false,
|
||||
nodes,
|
||||
alreadyRanIterations,
|
||||
hybrid;
|
||||
|
||||
if(simulation != undefined) initSimulation();
|
||||
if(forceS != undefined || forceF != undefined) initForces();
|
||||
|
||||
// Performed on first run
|
||||
function initialize() {
|
||||
initAlready = true;
|
||||
alreadyRanIterations = 0;
|
||||
simulation
|
||||
.on("tick", sampleTick)
|
||||
.on("end", sampleEnded)
|
||||
.nodes(sample)
|
||||
.force("Sample force", forceSample);
|
||||
console.log("Initialized Simulation for Hybrid");
|
||||
}
|
||||
|
||||
function initForces(){
|
||||
if (forceSample.onStableVelo) {
|
||||
forceSample.onStableVelo(sampleEnded);
|
||||
}
|
||||
|
||||
if (forceFull.onStableVelo) {
|
||||
forceFull.onStableVelo(fullEnded);
|
||||
}
|
||||
|
||||
// Set default value for interpDistanceFn if not been specified yet
|
||||
if(interpDistanceFn === undefined) {
|
||||
if(forceFull.distance == 'function')
|
||||
interpDistanceFn = forceFull.distance();
|
||||
else
|
||||
interpDistanceFn = constant(300);
|
||||
}
|
||||
}
|
||||
|
||||
function initSimulation(){
|
||||
nodes = simulation.nodes();
|
||||
simulation
|
||||
.stop()
|
||||
.alphaDecay(0)
|
||||
.alpha(1)
|
||||
|
||||
let sets = takeSampleFrom(nodes, Math.sqrt(nodes.length));
|
||||
sample = sets.sample;
|
||||
remainder = sets.remainder;
|
||||
}
|
||||
|
||||
// Sample simulation ticked 1 frame, keep track of number of iterations here.
|
||||
function sampleTick() {
|
||||
event.call("sampleTick");
|
||||
if(alreadyRanIterations++ >= SAMPLE_ITERATIONS){
|
||||
sampleEnded();
|
||||
}
|
||||
}
|
||||
|
||||
// Full simulation ticked 1 frame, keep track of number of iterations here.
|
||||
function fullTick() {
|
||||
event.call("fullTick");
|
||||
if(alreadyRanIterations++ >= FULL_ITERATIONS){
|
||||
fullEnded();
|
||||
}
|
||||
}
|
||||
|
||||
function fullEnded() {
|
||||
simulation.stop();
|
||||
initAlready = false;
|
||||
simulation.force("Full force", null);
|
||||
event.call("end");
|
||||
}
|
||||
|
||||
function sampleEnded() {
|
||||
simulation.stop();
|
||||
simulation.force("Sample force", null);
|
||||
// Reset velocity of all nodes
|
||||
for (let i=sample.length-1; i>=0; i--){
|
||||
sample[i].vx=0;
|
||||
sample[i].vy=0;
|
||||
}
|
||||
|
||||
event.call("startInterp");
|
||||
if (NUM_PIVOTS>=1) {
|
||||
interpolationPivots(sample, remainder, NUM_PIVOTS, interpDistanceFn, INTERP_FINE_ITS);
|
||||
} else {
|
||||
interpBruteForce(sample, remainder, interpDistanceFn, INTERP_FINE_ITS);
|
||||
}
|
||||
|
||||
event.call("fullTick");
|
||||
alreadyRanIterations = 0;
|
||||
simulation
|
||||
.on("tick", null)
|
||||
.on("end", null) // The ending condition should be iterations count
|
||||
.nodes(nodes);
|
||||
|
||||
if (FULL_ITERATIONS<1 || forceF === undefined || forceF === null) {
|
||||
event.call("end");
|
||||
return;
|
||||
}
|
||||
simulation
|
||||
.on("tick", fullTick)
|
||||
.force("Full force", forceFull)
|
||||
.restart();
|
||||
}
|
||||
|
||||
return hybrid = {
|
||||
restart: function () {
|
||||
if(!initAlready) initialize();
|
||||
simulation.restart();
|
||||
return hybrid;
|
||||
},
|
||||
|
||||
stop: function () {
|
||||
simulation.stop();
|
||||
return hybrid;
|
||||
},
|
||||
|
||||
numPivots: function (_) {
|
||||
return arguments.length ? (NUM_PIVOTS = +_, hybrid) : NUM_PIVOTS;
|
||||
},
|
||||
|
||||
sampleIterations: function (_) {
|
||||
return arguments.length ? (SAMPLE_ITERATIONS = +_, hybrid) : SAMPLE_ITERATIONS;
|
||||
},
|
||||
|
||||
fullIterations: function (_) {
|
||||
return arguments.length ? (FULL_ITERATIONS = +_, hybrid) : FULL_ITERATIONS;
|
||||
},
|
||||
|
||||
interpFindTuneIts: function (_) {
|
||||
return arguments.length ? (INTERP_FINE_ITS = +_, hybrid) : INTERP_FINE_ITS;
|
||||
},
|
||||
|
||||
on: function (name, _) {
|
||||
return arguments.length > 1 ? (event.on(name, _), hybrid) : event.on(name);
|
||||
},
|
||||
|
||||
subSet: function (_) {
|
||||
return arguments.length ? (sample = _, hybrid) : sample;
|
||||
},
|
||||
|
||||
nonSubSet: function (_) {
|
||||
return arguments.length ? (remainder = _, hybrid) : remainder;
|
||||
},
|
||||
|
||||
interpDistanceFn: function (_) {
|
||||
return arguments.length ? (interpDistanceFn = typeof _ === "function" ? _ : constant(+_), hybrid) : interpDistanceFn;
|
||||
},
|
||||
|
||||
simulation: function (_) {
|
||||
return arguments.length ? (initAlready = false, simulation = _, hybrid) : simulation;
|
||||
},
|
||||
|
||||
forceSample: function (_) {
|
||||
return arguments.length ? (forceSample = _, initForces(), hybrid) : forceSample;
|
||||
},
|
||||
|
||||
forceFull: function (_) {
|
||||
return arguments.length ? (forceFull = _, initForces(), hybrid) : forceFull;
|
||||
},
|
||||
|
||||
};
|
||||
}
|
||||
import {dispatch} from "d3-dispatch";
|
||||
import constant from "./constant";
|
||||
import interpBruteForce from "./interpolation/interpBruteForce";
|
||||
import interpolationPivots from "./interpolation/interpolationPivots";
|
||||
import {takeSampleFrom} from "./interpolation/helpers";
|
||||
|
||||
/**
|
||||
* An implementation of Chalmers, Morrison, and Ross' 2002 hybrid layout
|
||||
* algorithm with an option to use the 2003 pivot-based near neighbour searching
|
||||
* method.
|
||||
* It performs the 1996 neighbour sampling spring simulation model, on a
|
||||
* "sample set", sqrt(n) samples of the data.
|
||||
* Other data points are then interpolated into the model.
|
||||
* Finally, another spring simulation may be performed on the entire dataset to
|
||||
* clean up the model.
|
||||
* @param {object} sim - D3 Simulation object
|
||||
* @param {object} forceS - Pre-configured D3 force object for the sample set.
|
||||
The ending handler will be re-configured.
|
||||
Neighbour sampling force is expected, but other D3
|
||||
forces may also work.
|
||||
* @param {object} forceF - Pre-configured D3 force object for the simultion ran
|
||||
on the entire dataset at the end.
|
||||
Neighbour sampling force is expected, but other D3
|
||||
forces may also work.
|
||||
The force should not have any ending condition.
|
||||
*/
|
||||
export default function (sim, forceS, forceF) {
|
||||
var
|
||||
SAMPLE_ITERATIONS = 300,
|
||||
FULL_ITERATIONS = 20,
|
||||
interpDistanceFn,
|
||||
NUM_PIVOTS = 0,
|
||||
INTERP_FINE_ITS = 20,
|
||||
sample = [],
|
||||
remainder = [],
|
||||
simulation = sim,
|
||||
forceSample = forceS,
|
||||
forceFull = forceF,
|
||||
event = d3.dispatch("sampleTick", "fullTick", "startInterp", "end"),
|
||||
initAlready = false,
|
||||
nodes,
|
||||
alreadyRanIterations,
|
||||
hybrid;
|
||||
|
||||
if(simulation != undefined) initSimulation();
|
||||
if(forceS != undefined || forceF != undefined) initForces();
|
||||
|
||||
// Performed on first run
|
||||
function initialize() {
|
||||
initAlready = true;
|
||||
alreadyRanIterations = 0;
|
||||
simulation
|
||||
.on("tick", sampleTick)
|
||||
.on("end", sampleEnded)
|
||||
.nodes(sample)
|
||||
.force("Sample force", forceSample);
|
||||
console.log("Initialized Simulation for Hybrid");
|
||||
}
|
||||
|
||||
function initForces(){
|
||||
if (forceSample.onStableVelo) {
|
||||
forceSample.onStableVelo(sampleEnded);
|
||||
}
|
||||
|
||||
if (forceFull.onStableVelo) {
|
||||
forceFull.onStableVelo(fullEnded);
|
||||
}
|
||||
|
||||
// Set default value for interpDistanceFn if not been specified yet
|
||||
if(interpDistanceFn === undefined) {
|
||||
if(forceFull.distance == 'function')
|
||||
interpDistanceFn = forceFull.distance();
|
||||
else
|
||||
interpDistanceFn = constant(300);
|
||||
}
|
||||
}
|
||||
|
||||
function initSimulation(){
|
||||
nodes = simulation.nodes();
|
||||
simulation
|
||||
.stop()
|
||||
.alphaDecay(0)
|
||||
.alpha(1)
|
||||
|
||||
let sets = takeSampleFrom(nodes, Math.sqrt(nodes.length));
|
||||
sample = sets.sample;
|
||||
remainder = sets.remainder;
|
||||
}
|
||||
|
||||
// Sample simulation ticked 1 frame, keep track of number of iterations here.
|
||||
function sampleTick() {
|
||||
event.call("sampleTick");
|
||||
if(alreadyRanIterations++ >= SAMPLE_ITERATIONS){
|
||||
sampleEnded();
|
||||
}
|
||||
}
|
||||
|
||||
// Full simulation ticked 1 frame, keep track of number of iterations here.
|
||||
function fullTick() {
|
||||
event.call("fullTick");
|
||||
if(alreadyRanIterations++ >= FULL_ITERATIONS){
|
||||
fullEnded();
|
||||
}
|
||||
}
|
||||
|
||||
function fullEnded() {
|
||||
simulation.stop();
|
||||
initAlready = false;
|
||||
simulation.force("Full force", null);
|
||||
event.call("end");
|
||||
}
|
||||
|
||||
function sampleEnded() {
|
||||
simulation.stop();
|
||||
simulation.force("Sample force", null);
|
||||
// Reset velocity of all nodes
|
||||
for (let i=sample.length-1; i>=0; i--){
|
||||
sample[i].vx=0;
|
||||
sample[i].vy=0;
|
||||
}
|
||||
|
||||
event.call("startInterp");
|
||||
if (NUM_PIVOTS>=1) {
|
||||
interpolationPivots(sample, remainder, NUM_PIVOTS, interpDistanceFn, INTERP_FINE_ITS);
|
||||
} else {
|
||||
interpBruteForce(sample, remainder, interpDistanceFn, INTERP_FINE_ITS);
|
||||
}
|
||||
|
||||
event.call("fullTick");
|
||||
alreadyRanIterations = 0;
|
||||
simulation
|
||||
.on("tick", null)
|
||||
.on("end", null) // The ending condition should be iterations count
|
||||
.nodes(nodes);
|
||||
|
||||
if (FULL_ITERATIONS<1 || forceF === undefined || forceF === null) {
|
||||
event.call("end");
|
||||
return;
|
||||
}
|
||||
simulation
|
||||
.on("tick", fullTick)
|
||||
.force("Full force", forceFull)
|
||||
.restart();
|
||||
}
|
||||
|
||||
return hybrid = {
|
||||
restart: function () {
|
||||
if(!initAlready) initialize();
|
||||
simulation.restart();
|
||||
return hybrid;
|
||||
},
|
||||
|
||||
stop: function () {
|
||||
simulation.stop();
|
||||
return hybrid;
|
||||
},
|
||||
|
||||
numPivots: function (_) {
|
||||
return arguments.length ? (NUM_PIVOTS = +_, hybrid) : NUM_PIVOTS;
|
||||
},
|
||||
|
||||
sampleIterations: function (_) {
|
||||
return arguments.length ? (SAMPLE_ITERATIONS = +_, hybrid) : SAMPLE_ITERATIONS;
|
||||
},
|
||||
|
||||
fullIterations: function (_) {
|
||||
return arguments.length ? (FULL_ITERATIONS = +_, hybrid) : FULL_ITERATIONS;
|
||||
},
|
||||
|
||||
interpFindTuneIts: function (_) {
|
||||
return arguments.length ? (INTERP_FINE_ITS = +_, hybrid) : INTERP_FINE_ITS;
|
||||
},
|
||||
|
||||
on: function (name, _) {
|
||||
return arguments.length > 1 ? (event.on(name, _), hybrid) : event.on(name);
|
||||
},
|
||||
|
||||
subSet: function (_) {
|
||||
return arguments.length ? (sample = _, hybrid) : sample;
|
||||
},
|
||||
|
||||
nonSubSet: function (_) {
|
||||
return arguments.length ? (remainder = _, hybrid) : remainder;
|
||||
},
|
||||
|
||||
interpDistanceFn: function (_) {
|
||||
return arguments.length ? (interpDistanceFn = typeof _ === "function" ? _ : constant(+_), hybrid) : interpDistanceFn;
|
||||
},
|
||||
|
||||
simulation: function (_) {
|
||||
return arguments.length ? (initAlready = false, simulation = _, hybrid) : simulation;
|
||||
},
|
||||
|
||||
forceSample: function (_) {
|
||||
return arguments.length ? (forceSample = _, initForces(), hybrid) : forceSample;
|
||||
},
|
||||
|
||||
forceFull: function (_) {
|
||||
return arguments.length ? (forceFull = _, initForces(), hybrid) : forceFull;
|
||||
},
|
||||
|
||||
};
|
||||
}
|
||||
|
||||
@@ -10,8 +10,8 @@
|
||||
*/
|
||||
export function takeSampleFrom(sourceList, amount) {
|
||||
let randElements = [],
|
||||
max = sourceList.length,
|
||||
swap = false;
|
||||
max = sourceList.length,
|
||||
swap = false;
|
||||
|
||||
if (amount >= max) {
|
||||
return {sample: sourceList, remainder: {}};
|
||||
@@ -40,8 +40,7 @@ export function takeSampleFrom(sourceList, amount) {
|
||||
sample: remainder,
|
||||
remainder: randElements
|
||||
};
|
||||
}
|
||||
else {
|
||||
} else {
|
||||
return {
|
||||
sample: randElements,
|
||||
remainder: remainder
|
||||
|
||||
@@ -23,7 +23,7 @@ export default function(sampleSet, remainderSet, distanceFn, endingIts) {
|
||||
sampleSubset = takeSampleFrom(sampleSet, Math.sqrt(sampleSet.length)).sample,
|
||||
sampleSubsetDistanceCache = [];
|
||||
|
||||
// For each datapoint "node" to be interpolated
|
||||
// For each datapoint "node" to be interpolated
|
||||
for (let i = remainderSet.length-1; i>=0; i--) {
|
||||
let
|
||||
node = remainderSet[i],
|
||||
|
||||
@@ -36,7 +36,7 @@ export function placeNearToNearestNeighbour(node, nearNeighbour, radius, sampleS
|
||||
lowBound = 0.0,
|
||||
highBound = 0.0;
|
||||
|
||||
// Determine the closest quadrant
|
||||
// Determine the closest quadrant
|
||||
if (dist0 == dist180) {
|
||||
if (dist90 > dist270)
|
||||
lowBound = highBound = 270;
|
||||
@@ -55,14 +55,12 @@ export function placeNearToNearestNeighbour(node, nearNeighbour, radius, sampleS
|
||||
lowBound = 90;
|
||||
highBound = 180;
|
||||
}
|
||||
} else if (dist90 > dist270) {
|
||||
lowBound = 270;
|
||||
highBound = 360;
|
||||
} else {
|
||||
if (dist90 > dist270) {
|
||||
lowBound = 270;
|
||||
highBound = 360;
|
||||
} else {
|
||||
lowBound = 0;
|
||||
highBound = 90;
|
||||
}
|
||||
lowBound = 0;
|
||||
highBound = 90;
|
||||
}
|
||||
|
||||
// Determine the angle
|
||||
@@ -84,8 +82,8 @@ export function placeNearToNearestNeighbour(node, nearNeighbour, radius, sampleS
|
||||
|
||||
function sumForcesToSample(node, samples, sampleCache) {
|
||||
let nodeVx = 0,
|
||||
nodeVy = 0,
|
||||
x, y, l, i, sample;
|
||||
nodeVy = 0,
|
||||
x, y, l, i, sample;
|
||||
|
||||
for (i = samples.length-1; i >=0 ; i--) {
|
||||
sample = samples[i];
|
||||
|
||||
216
src/link.js
216
src/link.js
@@ -1,108 +1,108 @@
|
||||
import constant from "./constant";
|
||||
import jiggle from "./jiggle";
|
||||
|
||||
/**
|
||||
* Modified link force algorithm
|
||||
* - simplify calculations for parameters locked for spring model
|
||||
* - replace the use of links {} with loop. greatly reduce memory usage
|
||||
* - removed other unused functions
|
||||
* Alpha should be constant 1 for accurate simulation
|
||||
*/
|
||||
export default function() {
|
||||
var dataSizeFactor,
|
||||
distance = constant(30),
|
||||
distances = [],
|
||||
nodes,
|
||||
stableVelocity = 0,
|
||||
stableVeloHandler = null,
|
||||
latestVelocityDiff = 0,
|
||||
iterations = 1;
|
||||
|
||||
function force(alpha) {
|
||||
let n = nodes.length;
|
||||
// Cache old velocity for comparison later
|
||||
if (stableVeloHandler!==null && stableVelocity>=0) {
|
||||
for (let i = n-1, node; i>=0; i--) {
|
||||
node = nodes[i];
|
||||
node.oldvx = node.vx;
|
||||
node.oldvy = node.vy;
|
||||
}
|
||||
}
|
||||
|
||||
// Each iteration in a tick
|
||||
for (var k = 0, source, target, i, j, x, y, l; k < iterations; ++k) {
|
||||
// For each link
|
||||
for (i = 1; i < n; i++) for (j = 0; j < i; j++) {
|
||||
// jiggle so l won't be zero and divide by zero error after this
|
||||
source = nodes[i];
|
||||
target = nodes[j];
|
||||
x = target.x + target.vx - source.x - source.vx || jiggle();
|
||||
y = target.y + target.vy - source.y - source.vy || jiggle();
|
||||
l = Math.sqrt(x * x + y * y);
|
||||
l = (l - distances[i*(i-1)/2+j]) / l * dataSizeFactor * alpha;
|
||||
x *= l, y *= l;
|
||||
target.vx -= x;
|
||||
target.vy -= y;
|
||||
source.vx += x;
|
||||
source.vy += y;
|
||||
}
|
||||
}
|
||||
|
||||
// Calculate velocity changes, aka force applied.
|
||||
if (stableVeloHandler!==null && stableVelocity>=0) {
|
||||
let velocityDiff = 0;
|
||||
for (let i = n-1, node; i>=0; i--) {
|
||||
node = nodes[i];
|
||||
velocityDiff += Math.abs(Math.hypot(node.vx-node.oldvx, node.vy-node.oldvy));
|
||||
}
|
||||
velocityDiff /= n;
|
||||
latestVelocityDiff = velocityDiff;
|
||||
|
||||
if(velocityDiff<stableVelocity){
|
||||
stableVeloHandler();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
function initialize() {
|
||||
if (!nodes) return;
|
||||
// 0.5 to divide the force to two part for source and target node
|
||||
dataSizeFactor = 0.5/(nodes.length-1);
|
||||
initializeDistance();
|
||||
}
|
||||
|
||||
function initializeDistance() {
|
||||
if (!nodes) return;
|
||||
for (let i = 1, n = nodes.length; i < n; i++) {
|
||||
for (let j = 0; j < i; j++) {
|
||||
distances.push(distance(nodes[i], nodes[j]));
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
force.initialize = function(_) {
|
||||
nodes = _;
|
||||
initialize();
|
||||
};
|
||||
|
||||
force.iterations = function(_) {
|
||||
return arguments.length ? (iterations = +_, force) : iterations;
|
||||
};
|
||||
|
||||
force.distance = function(_) {
|
||||
return arguments.length ? (distance = typeof _ === "function" ? _ : constant(+_), initializeDistance(), force) : distance;
|
||||
};
|
||||
|
||||
force.latestAccel = function () {
|
||||
return latestVelocityDiff;
|
||||
};
|
||||
|
||||
force.onStableVelo = function (_) {
|
||||
return arguments.length ? (stableVeloHandler = _, force) : stableVeloHandler;
|
||||
};
|
||||
|
||||
force.stableVelocity = function (_) {
|
||||
return arguments.length ? (stableVelocity = _, force) : stableVelocity;
|
||||
};
|
||||
return force;
|
||||
}
|
||||
import constant from "./constant";
|
||||
import jiggle from "./jiggle";
|
||||
|
||||
/**
|
||||
* Modified link force algorithm
|
||||
* - simplify calculations for parameters locked for spring model
|
||||
* - replace the use of links {} with loop. greatly reduce memory usage
|
||||
* - removed other unused functions
|
||||
* Alpha should be constant 1 for accurate simulation
|
||||
*/
|
||||
export default function() {
|
||||
var dataSizeFactor,
|
||||
distance = constant(30),
|
||||
distances = [],
|
||||
nodes,
|
||||
stableVelocity = 0,
|
||||
stableVeloHandler = null,
|
||||
latestVelocityDiff = 0,
|
||||
iterations = 1;
|
||||
|
||||
function force(alpha) {
|
||||
let n = nodes.length;
|
||||
// Cache old velocity for comparison later
|
||||
if (stableVeloHandler!==null && stableVelocity>=0) {
|
||||
for (let i = n-1, node; i>=0; i--) {
|
||||
node = nodes[i];
|
||||
node.oldvx = node.vx;
|
||||
node.oldvy = node.vy;
|
||||
}
|
||||
}
|
||||
|
||||
// Each iteration in a tick
|
||||
for (var k = 0, source, target, i, j, x, y, l; k < iterations; ++k) {
|
||||
// For each link
|
||||
for (i = 1; i < n; i++) for (j = 0; j < i; j++) {
|
||||
// jiggle so l won't be zero and divide by zero error after this
|
||||
source = nodes[i];
|
||||
target = nodes[j];
|
||||
x = target.x + target.vx - source.x - source.vx || jiggle();
|
||||
y = target.y + target.vy - source.y - source.vy || jiggle();
|
||||
l = Math.sqrt(x * x + y * y);
|
||||
l = (l - distances[i*(i-1)/2+j]) / l * dataSizeFactor * alpha;
|
||||
x *= l, y *= l;
|
||||
target.vx -= x;
|
||||
target.vy -= y;
|
||||
source.vx += x;
|
||||
source.vy += y;
|
||||
}
|
||||
}
|
||||
|
||||
// Calculate velocity changes, aka force applied.
|
||||
if (stableVeloHandler!==null && stableVelocity>=0) {
|
||||
let velocityDiff = 0;
|
||||
for (let i = n-1, node; i>=0; i--) {
|
||||
node = nodes[i];
|
||||
velocityDiff += Math.abs(Math.hypot(node.vx-node.oldvx, node.vy-node.oldvy));
|
||||
}
|
||||
velocityDiff /= n;
|
||||
latestVelocityDiff = velocityDiff;
|
||||
|
||||
if(velocityDiff<stableVelocity){
|
||||
stableVeloHandler();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
function initialize() {
|
||||
if (!nodes) return;
|
||||
// 0.5 to divide the force to two part for source and target node
|
||||
dataSizeFactor = 0.5/(nodes.length-1);
|
||||
initializeDistance();
|
||||
}
|
||||
|
||||
function initializeDistance() {
|
||||
if (!nodes) return;
|
||||
for (let i = 1, n = nodes.length; i < n; i++) {
|
||||
for (let j = 0; j < i; j++) {
|
||||
distances.push(distance(nodes[i], nodes[j]));
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
force.initialize = function(_) {
|
||||
nodes = _;
|
||||
initialize();
|
||||
};
|
||||
|
||||
force.iterations = function(_) {
|
||||
return arguments.length ? (iterations = +_, force) : iterations;
|
||||
};
|
||||
|
||||
force.distance = function(_) {
|
||||
return arguments.length ? (distance = typeof _ === "function" ? _ : constant(+_), initializeDistance(), force) : distance;
|
||||
};
|
||||
|
||||
force.latestAccel = function () {
|
||||
return latestVelocityDiff;
|
||||
};
|
||||
|
||||
force.onStableVelo = function (_) {
|
||||
return arguments.length ? (stableVeloHandler = _, force) : stableVeloHandler;
|
||||
};
|
||||
|
||||
force.stableVelocity = function (_) {
|
||||
return arguments.length ? (stableVelocity = _, force) : stableVelocity;
|
||||
};
|
||||
return force;
|
||||
}
|
||||
|
||||
@@ -21,7 +21,7 @@ export default function () {
|
||||
dataSizeFactor,
|
||||
latestVelocityDiff = 0;
|
||||
|
||||
/**
|
||||
/**
|
||||
* Apply spring forces at each simulation iteration.
|
||||
* @param {number} alpha - multiplier for amount of force applied
|
||||
*/
|
||||
|
||||
752
src/t-sne.js
752
src/t-sne.js
@@ -1,375 +1,377 @@
|
||||
/* eslint-disable block-scoped-var */
|
||||
import constant from "./constant";
|
||||
|
||||
/**
|
||||
* Set the node id accessor to the specified i.
|
||||
* @param {node} d - node.
|
||||
* @param {accessor} i - id accessor.
|
||||
* @return {accessor} - node id accessor.
|
||||
*/
|
||||
function index(d, i) {
|
||||
return i;
|
||||
}
|
||||
|
||||
/**
|
||||
* t-SNE implementation in D3 by using the code existing in tsnejs
|
||||
* (https://github.com/karpathy/tsnejs) to compute the solution.
|
||||
*/
|
||||
export default function() {
|
||||
var id = index,
|
||||
distance = constant(300),
|
||||
nodes,
|
||||
perplexity = 30,
|
||||
learningRate = 10,
|
||||
iteration = 0,
|
||||
dim = 2,
|
||||
N, // length of the nodes.
|
||||
P, // probability matrix.
|
||||
Y, // solution.
|
||||
gains,
|
||||
ystep;
|
||||
|
||||
/**
|
||||
* Make a step in t-SNE algorithm and set the velocities for the nodes
|
||||
* to accumulate the values from solution.
|
||||
*/
|
||||
function force() {
|
||||
// Make a step at each iteration.
|
||||
step();
|
||||
var solution = getSolution();
|
||||
|
||||
// Set the velocity for each node using the solution.
|
||||
for (var i = 0; i < nodes.length; i++) {
|
||||
nodes[i].vx += solution[i][0];
|
||||
nodes[i].vy += solution[i][1];
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculates the random number from Gaussian distribution.
|
||||
* @return {number} random number.
|
||||
*/
|
||||
function gaussRandom() {
|
||||
let u = 2 * Math.random() - 1;
|
||||
let v = 2 * Math.random() - 1;
|
||||
let r = u * u + v * v;
|
||||
if (r == 0 || r > 1) {return gaussRandom();}
|
||||
return u * Math.sqrt(-2 * Math.log(r) / r);
|
||||
}
|
||||
|
||||
/**
|
||||
* Return the normalized number.
|
||||
* @return {number} normalized random number from Gaussian distribution.
|
||||
*/
|
||||
function randomN() {
|
||||
return gaussRandom() * 1e-4;
|
||||
}
|
||||
|
||||
function sign(x) {
|
||||
return x > 0 ? 1 : x < 0 ? -1 : 0;
|
||||
}
|
||||
|
||||
/**
|
||||
* Create an array of length n filled with zeros.
|
||||
* @param {number} n - length of array.
|
||||
* @return {Float64Array} - array of zeros with length n.
|
||||
*/
|
||||
function zeros(n) {
|
||||
if (typeof n === 'undefined' || isNaN(n)) {
|
||||
return [];
|
||||
}
|
||||
return new Float64Array(n); // typed arrays are faster
|
||||
}
|
||||
|
||||
// Returns a 2d array of random numbers
|
||||
/**
|
||||
* Creates a 2d array filled with random numbers.
|
||||
* @param {number} n - rows.
|
||||
* @param {number} d - columns.
|
||||
* @return {array} - 2d array
|
||||
*/
|
||||
function random2d(n, d) {
|
||||
var x = [];
|
||||
for (var i = 0; i < n; i++) {
|
||||
var y = [];
|
||||
for (var j = 0; j < d; j++) {
|
||||
y.push(randomN());
|
||||
}
|
||||
x.push(y);
|
||||
}
|
||||
return x;
|
||||
}
|
||||
|
||||
/**
|
||||
* Compute the probability matrix using the provided data.
|
||||
* @param {array} data - nodes.
|
||||
* @param {number} perplexity - used to calculate entropy of distribution.
|
||||
* @param {number} tol - limit for entropy difference.
|
||||
* @return {2d array} - 2d matrix containing probabilities.
|
||||
*/
|
||||
function d2p(data, perplexity, tol) {
|
||||
N = Math.floor(data.length);
|
||||
var Htarget = Math.log(perplexity); // target entropy of distribution.
|
||||
var P1 = zeros(N * N); // temporary probability matrix.
|
||||
|
||||
var prow = zeros(N); // a temporary storage compartment.
|
||||
for (var i = 0; i < N; i++) {
|
||||
var betamin = -Infinity;
|
||||
var betamax = Infinity;
|
||||
var beta = 1; // initial value of precision.
|
||||
var done = false;
|
||||
var maxtries = 50;
|
||||
|
||||
// Perform binary search to find a suitable precision beta
|
||||
// so that the entropy of the distribution is appropriate.
|
||||
var num = 0;
|
||||
while (!done) {
|
||||
// Compute entropy and kernel row with beta precision.
|
||||
var psum = 0.0;
|
||||
for (var j = 0; j < N; j++) {
|
||||
var pj = Math.exp(-distance(data[i], data[j]) * beta);
|
||||
// Ignore the diagonals
|
||||
if (i === j) {
|
||||
pj = 0;
|
||||
}
|
||||
prow[j] = pj;
|
||||
psum += pj;
|
||||
}
|
||||
// Normalize p and compute entropy.
|
||||
var Hhere = 0.0;
|
||||
for (j = 0; j < N; j++) {
|
||||
if (psum == 0) {
|
||||
pj = 0;
|
||||
} else {
|
||||
pj = prow[j] / psum;
|
||||
}
|
||||
prow[j] = pj;
|
||||
if (pj > 1e-7) {
|
||||
Hhere -= pj * Math.log(pj);
|
||||
}
|
||||
}
|
||||
|
||||
// Adjust beta based on result.
|
||||
if (Hhere > Htarget) {
|
||||
// Entropy was too high (distribution too diffuse)
|
||||
// so we need to increase the precision for more peaky distribution.
|
||||
betamin = beta; // move up the bounds.
|
||||
if (betamax === Infinity) {
|
||||
beta = beta * 2;
|
||||
} else {
|
||||
beta = (beta + betamax) / 2;
|
||||
}
|
||||
|
||||
} else {
|
||||
// Converse case. Make distrubtion less peaky.
|
||||
betamax = beta;
|
||||
if (betamin === -Infinity) {
|
||||
beta = beta / 2;
|
||||
} else {
|
||||
beta = (beta + betamin) / 2;
|
||||
}
|
||||
}
|
||||
|
||||
// Stopping conditions: too many tries or got a good precision.
|
||||
num++;
|
||||
if (Math.abs(Hhere - Htarget) < tol || num >= maxtries) {
|
||||
done = true;
|
||||
}
|
||||
}
|
||||
|
||||
// Copy over the final prow to P1 at row i
|
||||
for (j = 0; j < N; j++) {
|
||||
P1[i * N + j] = prow[j];
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
// Symmetrize P and normalize it to sum to 1 over all ij
|
||||
var Pout = zeros(N * N);
|
||||
var N2 = N * 2;
|
||||
for (i = 0; i < N; i++) {
|
||||
for (j = 0; j < N; j++) {
|
||||
Pout[i * N + j] = Math.max((P1[i * N + j] + P1[j * N + i]) / N2, 1e-100);
|
||||
}
|
||||
}
|
||||
return Pout;
|
||||
}
|
||||
|
||||
/**
|
||||
* Initialize a starting (random) solution.
|
||||
*/
|
||||
function initSolution() {
|
||||
Y = random2d(N, dim);
|
||||
// Step gains to accelerate progress in unchanging directions.
|
||||
gains = random2d(N, dim, 1.0);
|
||||
// Momentum accumulator.
|
||||
ystep = random2d(N, dim, 0.0);
|
||||
iteration = 0;
|
||||
}
|
||||
|
||||
/**
|
||||
* @return {2d array} the solution.
|
||||
*/
|
||||
function getSolution() {
|
||||
return Y;
|
||||
}
|
||||
|
||||
/**
|
||||
* Do a single step (iteration) for the layout.
|
||||
* @return {number} the current cost.
|
||||
*/
|
||||
function step() {
|
||||
iteration += 1;
|
||||
|
||||
var cg = costGrad(Y); // Evaluate gradient.
|
||||
var cost = cg.cost;
|
||||
var grad = cg.grad;
|
||||
|
||||
// Perform gradient step.
|
||||
var ymean = zeros(dim);
|
||||
for (var i = 0; i < N; i++) {
|
||||
for (var d = 0; d < dim; d++) {
|
||||
var gid = grad[i][d];
|
||||
var sid = ystep[i][d];
|
||||
var gainid = gains[i][d];
|
||||
|
||||
// Compute gain update.
|
||||
var newgain = sign(gid) === sign(sid) ? gainid * 0.8 : gainid + 0.2;
|
||||
if (newgain < 0.01) {
|
||||
newgain = 0.01;
|
||||
}
|
||||
gains[i][d] = newgain;
|
||||
|
||||
// Compute momentum step direction.
|
||||
var momval = iteration < 250 ? 0.5 : 0.8;
|
||||
var newsid = momval * sid - learningRate * newgain * grad[i][d];
|
||||
ystep[i][d] = newsid;
|
||||
|
||||
// Do the step.
|
||||
Y[i][d] += newsid;
|
||||
|
||||
// Accumulate mean so that we can center later.
|
||||
ymean[d] += Y[i][d];
|
||||
}
|
||||
}
|
||||
|
||||
// Reproject Y to have the zero mean.
|
||||
for (i = 0; i < N; i++) {
|
||||
for (d = 0; d < dim; d++) {
|
||||
Y[i][d] -= ymean[d] / N;
|
||||
}
|
||||
}
|
||||
return cost;
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculate the cost and the gradient.
|
||||
* @param {2d array} Y - the current solution to evaluate.
|
||||
* @return {object} that contains a cost and a gradient.
|
||||
*/
|
||||
function costGrad(Y) {
|
||||
|
||||
var pmul = iteration < 100 ? 4 : 1;
|
||||
|
||||
// Compute current Q distribution, unnormalized first.
|
||||
var Qu = zeros(N * N);
|
||||
var qsum = 0.0;
|
||||
for (var i = 0; i < N; i++) {
|
||||
for (var j = i + 1; j < N; j++) {
|
||||
var dsum = 0.0;
|
||||
for (var d = 0; d < dim; d++) {
|
||||
var dhere = Y[i][d] - Y[j][d];
|
||||
dsum += dhere * dhere;
|
||||
}
|
||||
var qu = 1.0 / (1.0 + dsum); // Student t-distribution.
|
||||
Qu[i * N + j] = qu;
|
||||
Qu[j * N + i] = qu;
|
||||
qsum += 2 * qu;
|
||||
}
|
||||
}
|
||||
// Normalize Q distribution to sum to 1.
|
||||
var NN = N * N;
|
||||
var Q = zeros(NN);
|
||||
for (var q = 0; q < NN; q++) {
|
||||
Q[q] = Math.max(Qu[q] / qsum, 1e-100);
|
||||
}
|
||||
|
||||
var cost = 0.0;
|
||||
var grad = [];
|
||||
for (i = 0; i < N; i++) {
|
||||
var gsum = new Array(dim); // Initialize gradiet for point i.
|
||||
for (d = 0; d < dim; d++) {
|
||||
gsum[d] = 0.0;
|
||||
}
|
||||
for (j = 0; j < N; j++) {
|
||||
// Accumulate the cost.
|
||||
cost += -P[i * N + j] * Math.log(Q[i * N + j]);
|
||||
var premult = 4 * (pmul * P[i * N + j] - Q[i * N + j]) * Qu[i * N + j];
|
||||
for (d = 0; d < dim; d++) {
|
||||
gsum[d] += premult * (Y[i][d] - Y[j][d]);
|
||||
}
|
||||
}
|
||||
grad.push(gsum);
|
||||
}
|
||||
|
||||
return {
|
||||
cost: cost,
|
||||
grad: grad
|
||||
};
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculates the stress. Basically, it computes the difference between
|
||||
* high dimensional distance and real distance. The lower the stress is,
|
||||
* the better layout.
|
||||
* @return {number} - stress of the layout.
|
||||
*/
|
||||
function getStress() {
|
||||
var totalDiffSq = 0,
|
||||
totalHighDistSq = 0;
|
||||
for (var i = 0, source, target, realDist, highDist; i < nodes.length; i++) {
|
||||
for (var j = 0; j < nodes.length; j++) {
|
||||
if (i !== j) {
|
||||
source = nodes[i], target = nodes[j];
|
||||
realDist = Math.hypot(target.x - source.x, target.y - source.y);
|
||||
highDist = +distance(nodes[i], nodes[j]);
|
||||
totalDiffSq += Math.pow(realDist - highDist, 2);
|
||||
totalHighDistSq += highDist * highDist;
|
||||
}
|
||||
}
|
||||
}
|
||||
return Math.sqrt(totalDiffSq / totalHighDistSq);
|
||||
}
|
||||
|
||||
// API for initializing the algorithm, setting parameters and querying
|
||||
// metrics.
|
||||
force.initialize = function(_) {
|
||||
nodes = _;
|
||||
N = nodes.length;
|
||||
// Initialize the probability matrix.
|
||||
P = d2p(nodes, perplexity, 1e-4);
|
||||
initSolution();
|
||||
};
|
||||
|
||||
force.id = function(_) {
|
||||
return arguments.length ? (id = _, force) : id;
|
||||
};
|
||||
|
||||
force.distance = function(_) {
|
||||
return arguments.length ? (distance = typeof _ === "function" ? _ : constant(+_), force) : distance;
|
||||
};
|
||||
|
||||
force.stress = function() {
|
||||
return getStress();
|
||||
};
|
||||
|
||||
force.learningRate = function(_) {
|
||||
return arguments.length ? (learningRate = +_, force) : learningRate;
|
||||
};
|
||||
|
||||
force.perplexity = function(_) {
|
||||
return arguments.length ? (perplexity = +_, force) : perplexity;
|
||||
};
|
||||
|
||||
return force;
|
||||
}
|
||||
/* eslint-disable block-scoped-var */
|
||||
import constant from "./constant";
|
||||
|
||||
/**
|
||||
* Set the node id accessor to the specified i.
|
||||
* @param {node} d - node.
|
||||
* @param {accessor} i - id accessor.
|
||||
* @return {accessor} - node id accessor.
|
||||
*/
|
||||
function index(d, i) {
|
||||
return i;
|
||||
}
|
||||
|
||||
/**
|
||||
* t-SNE implementation in D3 by using the code existing in tsnejs
|
||||
* (https://github.com/karpathy/tsnejs) to compute the solution.
|
||||
*/
|
||||
export default function() {
|
||||
var id = index,
|
||||
distance = constant(300),
|
||||
nodes,
|
||||
perplexity = 30,
|
||||
learningRate = 10,
|
||||
iteration = 0,
|
||||
dim = 2,
|
||||
N, // length of the nodes.
|
||||
P, // probability matrix.
|
||||
Y, // solution.
|
||||
gains,
|
||||
ystep;
|
||||
|
||||
/**
|
||||
* Make a step in t-SNE algorithm and set the velocities for the nodes
|
||||
* to accumulate the values from solution.
|
||||
*/
|
||||
function force() {
|
||||
// Make a step at each iteration.
|
||||
step();
|
||||
var solution = getSolution();
|
||||
|
||||
// Set the velocity for each node using the solution.
|
||||
for (var i = 0; i < nodes.length; i++) {
|
||||
nodes[i].vx += solution[i][0];
|
||||
nodes[i].vy += solution[i][1];
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculates the random number from Gaussian distribution.
|
||||
* @return {number} random number.
|
||||
*/
|
||||
function gaussRandom() {
|
||||
let u = 2 * Math.random() - 1;
|
||||
let v = 2 * Math.random() - 1;
|
||||
let r = u * u + v * v;
|
||||
if (r == 0 || r > 1) {
|
||||
return gaussRandom();
|
||||
}
|
||||
return u * Math.sqrt(-2 * Math.log(r) / r);
|
||||
}
|
||||
|
||||
/**
|
||||
* Return the normalized number.
|
||||
* @return {number} normalized random number from Gaussian distribution.
|
||||
*/
|
||||
function randomN() {
|
||||
return gaussRandom() * 1e-4;
|
||||
}
|
||||
|
||||
function sign(x) {
|
||||
return x > 0 ? 1 : x < 0 ? -1 : 0;
|
||||
}
|
||||
|
||||
/**
|
||||
* Create an array of length n filled with zeros.
|
||||
* @param {number} n - length of array.
|
||||
* @return {Float64Array} - array of zeros with length n.
|
||||
*/
|
||||
function zeros(n) {
|
||||
if (typeof n === 'undefined' || isNaN(n)) {
|
||||
return [];
|
||||
}
|
||||
return new Float64Array(n); // typed arrays are faster
|
||||
}
|
||||
|
||||
// Returns a 2d array of random numbers
|
||||
/**
|
||||
* Creates a 2d array filled with random numbers.
|
||||
* @param {number} n - rows.
|
||||
* @param {number} d - columns.
|
||||
* @return {array} - 2d array
|
||||
*/
|
||||
function random2d(n, d) {
|
||||
var x = [];
|
||||
for (var i = 0; i < n; i++) {
|
||||
var y = [];
|
||||
for (var j = 0; j < d; j++) {
|
||||
y.push(randomN());
|
||||
}
|
||||
x.push(y);
|
||||
}
|
||||
return x;
|
||||
}
|
||||
|
||||
/**
|
||||
* Compute the probability matrix using the provided data.
|
||||
* @param {array} data - nodes.
|
||||
* @param {number} perplexity - used to calculate entropy of distribution.
|
||||
* @param {number} tol - limit for entropy difference.
|
||||
* @return {2d array} - 2d matrix containing probabilities.
|
||||
*/
|
||||
function d2p(data, perplexity, tol) {
|
||||
N = Math.floor(data.length);
|
||||
var Htarget = Math.log(perplexity); // target entropy of distribution.
|
||||
var P1 = zeros(N * N); // temporary probability matrix.
|
||||
|
||||
var prow = zeros(N); // a temporary storage compartment.
|
||||
for (var i = 0; i < N; i++) {
|
||||
var betamin = -Infinity;
|
||||
var betamax = Infinity;
|
||||
var beta = 1; // initial value of precision.
|
||||
var done = false;
|
||||
var maxtries = 50;
|
||||
|
||||
// Perform binary search to find a suitable precision beta
|
||||
// so that the entropy of the distribution is appropriate.
|
||||
var num = 0;
|
||||
while (!done) {
|
||||
// Compute entropy and kernel row with beta precision.
|
||||
var psum = 0.0;
|
||||
for (var j = 0; j < N; j++) {
|
||||
var pj = Math.exp(-distance(data[i], data[j]) * beta);
|
||||
// Ignore the diagonals
|
||||
if (i === j) {
|
||||
pj = 0;
|
||||
}
|
||||
prow[j] = pj;
|
||||
psum += pj;
|
||||
}
|
||||
// Normalize p and compute entropy.
|
||||
var Hhere = 0.0;
|
||||
for (j = 0; j < N; j++) {
|
||||
if (psum == 0) {
|
||||
pj = 0;
|
||||
} else {
|
||||
pj = prow[j] / psum;
|
||||
}
|
||||
prow[j] = pj;
|
||||
if (pj > 1e-7) {
|
||||
Hhere -= pj * Math.log(pj);
|
||||
}
|
||||
}
|
||||
|
||||
// Adjust beta based on result.
|
||||
if (Hhere > Htarget) {
|
||||
// Entropy was too high (distribution too diffuse)
|
||||
// so we need to increase the precision for more peaky distribution.
|
||||
betamin = beta; // move up the bounds.
|
||||
if (betamax === Infinity) {
|
||||
beta = beta * 2;
|
||||
} else {
|
||||
beta = (beta + betamax) / 2;
|
||||
}
|
||||
|
||||
} else {
|
||||
// Converse case. Make distrubtion less peaky.
|
||||
betamax = beta;
|
||||
if (betamin === -Infinity) {
|
||||
beta = beta / 2;
|
||||
} else {
|
||||
beta = (beta + betamin) / 2;
|
||||
}
|
||||
}
|
||||
|
||||
// Stopping conditions: too many tries or got a good precision.
|
||||
num++;
|
||||
if (Math.abs(Hhere - Htarget) < tol || num >= maxtries) {
|
||||
done = true;
|
||||
}
|
||||
}
|
||||
|
||||
// Copy over the final prow to P1 at row i
|
||||
for (j = 0; j < N; j++) {
|
||||
P1[i * N + j] = prow[j];
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
// Symmetrize P and normalize it to sum to 1 over all ij
|
||||
var Pout = zeros(N * N);
|
||||
var N2 = N * 2;
|
||||
for (i = 0; i < N; i++) {
|
||||
for (j = 0; j < N; j++) {
|
||||
Pout[i * N + j] = Math.max((P1[i * N + j] + P1[j * N + i]) / N2, 1e-100);
|
||||
}
|
||||
}
|
||||
return Pout;
|
||||
}
|
||||
|
||||
/**
|
||||
* Initialize a starting (random) solution.
|
||||
*/
|
||||
function initSolution() {
|
||||
Y = random2d(N, dim);
|
||||
// Step gains to accelerate progress in unchanging directions.
|
||||
gains = random2d(N, dim, 1.0);
|
||||
// Momentum accumulator.
|
||||
ystep = random2d(N, dim, 0.0);
|
||||
iteration = 0;
|
||||
}
|
||||
|
||||
/**
|
||||
* @return {2d array} the solution.
|
||||
*/
|
||||
function getSolution() {
|
||||
return Y;
|
||||
}
|
||||
|
||||
/**
|
||||
* Do a single step (iteration) for the layout.
|
||||
* @return {number} the current cost.
|
||||
*/
|
||||
function step() {
|
||||
iteration += 1;
|
||||
|
||||
var cg = costGrad(Y); // Evaluate gradient.
|
||||
var cost = cg.cost;
|
||||
var grad = cg.grad;
|
||||
|
||||
// Perform gradient step.
|
||||
var ymean = zeros(dim);
|
||||
for (var i = 0; i < N; i++) {
|
||||
for (var d = 0; d < dim; d++) {
|
||||
var gid = grad[i][d];
|
||||
var sid = ystep[i][d];
|
||||
var gainid = gains[i][d];
|
||||
|
||||
// Compute gain update.
|
||||
var newgain = sign(gid) === sign(sid) ? gainid * 0.8 : gainid + 0.2;
|
||||
if (newgain < 0.01) {
|
||||
newgain = 0.01;
|
||||
}
|
||||
gains[i][d] = newgain;
|
||||
|
||||
// Compute momentum step direction.
|
||||
var momval = iteration < 250 ? 0.5 : 0.8;
|
||||
var newsid = momval * sid - learningRate * newgain * grad[i][d];
|
||||
ystep[i][d] = newsid;
|
||||
|
||||
// Do the step.
|
||||
Y[i][d] += newsid;
|
||||
|
||||
// Accumulate mean so that we can center later.
|
||||
ymean[d] += Y[i][d];
|
||||
}
|
||||
}
|
||||
|
||||
// Reproject Y to have the zero mean.
|
||||
for (i = 0; i < N; i++) {
|
||||
for (d = 0; d < dim; d++) {
|
||||
Y[i][d] -= ymean[d] / N;
|
||||
}
|
||||
}
|
||||
return cost;
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculate the cost and the gradient.
|
||||
* @param {2d array} Y - the current solution to evaluate.
|
||||
* @return {object} that contains a cost and a gradient.
|
||||
*/
|
||||
function costGrad(Y) {
|
||||
|
||||
var pmul = iteration < 100 ? 4 : 1;
|
||||
|
||||
// Compute current Q distribution, unnormalized first.
|
||||
var Qu = zeros(N * N);
|
||||
var qsum = 0.0;
|
||||
for (var i = 0; i < N; i++) {
|
||||
for (var j = i + 1; j < N; j++) {
|
||||
var dsum = 0.0;
|
||||
for (var d = 0; d < dim; d++) {
|
||||
var dhere = Y[i][d] - Y[j][d];
|
||||
dsum += dhere * dhere;
|
||||
}
|
||||
var qu = 1.0 / (1.0 + dsum); // Student t-distribution.
|
||||
Qu[i * N + j] = qu;
|
||||
Qu[j * N + i] = qu;
|
||||
qsum += 2 * qu;
|
||||
}
|
||||
}
|
||||
// Normalize Q distribution to sum to 1.
|
||||
var NN = N * N;
|
||||
var Q = zeros(NN);
|
||||
for (var q = 0; q < NN; q++) {
|
||||
Q[q] = Math.max(Qu[q] / qsum, 1e-100);
|
||||
}
|
||||
|
||||
var cost = 0.0;
|
||||
var grad = [];
|
||||
for (i = 0; i < N; i++) {
|
||||
var gsum = new Array(dim); // Initialize gradiet for point i.
|
||||
for (d = 0; d < dim; d++) {
|
||||
gsum[d] = 0.0;
|
||||
}
|
||||
for (j = 0; j < N; j++) {
|
||||
// Accumulate the cost.
|
||||
cost += -P[i * N + j] * Math.log(Q[i * N + j]);
|
||||
var premult = 4 * (pmul * P[i * N + j] - Q[i * N + j]) * Qu[i * N + j];
|
||||
for (d = 0; d < dim; d++) {
|
||||
gsum[d] += premult * (Y[i][d] - Y[j][d]);
|
||||
}
|
||||
}
|
||||
grad.push(gsum);
|
||||
}
|
||||
|
||||
return {
|
||||
cost: cost,
|
||||
grad: grad
|
||||
};
|
||||
}
|
||||
|
||||
/**
|
||||
* Calculates the stress. Basically, it computes the difference between
|
||||
* high dimensional distance and real distance. The lower the stress is,
|
||||
* the better layout.
|
||||
* @return {number} - stress of the layout.
|
||||
*/
|
||||
function getStress() {
|
||||
var totalDiffSq = 0,
|
||||
totalHighDistSq = 0;
|
||||
for (var i = 0, source, target, realDist, highDist; i < nodes.length; i++) {
|
||||
for (var j = 0; j < nodes.length; j++) {
|
||||
if (i !== j) {
|
||||
source = nodes[i], target = nodes[j];
|
||||
realDist = Math.hypot(target.x - source.x, target.y - source.y);
|
||||
highDist = +distance(nodes[i], nodes[j]);
|
||||
totalDiffSq += Math.pow(realDist - highDist, 2);
|
||||
totalHighDistSq += highDist * highDist;
|
||||
}
|
||||
}
|
||||
}
|
||||
return Math.sqrt(totalDiffSq / totalHighDistSq);
|
||||
}
|
||||
|
||||
// API for initializing the algorithm, setting parameters and querying
|
||||
// metrics.
|
||||
force.initialize = function(_) {
|
||||
nodes = _;
|
||||
N = nodes.length;
|
||||
// Initialize the probability matrix.
|
||||
P = d2p(nodes, perplexity, 1e-4);
|
||||
initSolution();
|
||||
};
|
||||
|
||||
force.id = function(_) {
|
||||
return arguments.length ? (id = _, force) : id;
|
||||
};
|
||||
|
||||
force.distance = function(_) {
|
||||
return arguments.length ? (distance = typeof _ === "function" ? _ : constant(+_), force) : distance;
|
||||
};
|
||||
|
||||
force.stress = function() {
|
||||
return getStress();
|
||||
};
|
||||
|
||||
force.learningRate = function(_) {
|
||||
return arguments.length ? (learningRate = +_, force) : learningRate;
|
||||
};
|
||||
|
||||
force.perplexity = function(_) {
|
||||
return arguments.length ? (perplexity = +_, force) : perplexity;
|
||||
};
|
||||
|
||||
return force;
|
||||
}
|
||||
|
||||
26
src/xy.js
26
src/xy.js
@@ -1,13 +1,13 @@
|
||||
/**
|
||||
* @return x value of a node
|
||||
*/
|
||||
export function x(d) {
|
||||
return d.x;
|
||||
}
|
||||
|
||||
/**
|
||||
* @return y value of a node
|
||||
*/
|
||||
export function y(d) {
|
||||
return d.y;
|
||||
}
|
||||
/**
|
||||
* @return x value of a node
|
||||
*/
|
||||
export function x(d) {
|
||||
return d.x;
|
||||
}
|
||||
|
||||
/**
|
||||
* @return y value of a node
|
||||
*/
|
||||
export function y(d) {
|
||||
return d.y;
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user