Published
Edited
May 21, 2020
13 stars
Jacobi Elliptic Functions in the Complex PlaneComputing π with the Bailey-Borwein-Plouffe FormulaAdaptive Contouring in Fragment ShadersDomain Coloring with Adaptive ContouringDouble Compound Pendulums3D Reaction-DiffusionQuasicrystals2D (Non-physical) N-body Gravity with Poisson's EquationPeriodic Planar Three-Body OrbitseasesFinding Roots in the Complex Planeglsl-complexglsl-complex-auto-differentiationpolynomial-rootscomplex-auto-differentiationutilsIntegration
Grouping Points with Principal Component Analysis
Also listed in…
More or less finished notebooks
axes = Eigen.ensure(invalidation, () => {
// `ensure` is a small memory management wrapper that cleans up resource
// allocated on the wasm heap (e.g. Matrixf in the return values) when
// this cell is invalidated.
// Accumulate the centroid:
let xc = 0;
let yc = 0;
// Components of cross-covariance matrix A
let A00 = 0;
let A10 = 0;
let A11 = 0;
// Iterate over points
let xi, yi;
for (let i = 0; i < n; i++) {
xi = x.get(i, 0);
yi = x.get(i, 1);
// Accumulate the centroid
xc += xi;
yc += yi;
// Sum moments
// NOTE: Don't do this! You should compute the centroid in one pass,
// then sum, e.g. (xi - xc) * (yi - yc) *relative to the centroid* in
// order to prevent catastrophic cancellation below.
A00 += xi * xi;
A10 += xi * yi;
A11 += yi * yi;
}
// Divide to get the centroid
xc /= n;
yc /= n;
// Unshift the cross-covariance by the centroid using the parallel axis theorem
// This may cause unnecessary catastrophic cancellation if the centroid is not
// right at the origin!
A00 = A00 / n - xc * xc;
A10 = A10 / n - yc * xc;
A11 = A11 / n - yc * yc;
return {
centroid: Matrixf.fromArray([xc, yc], [1, 2]),
covariance: Matrixf.fromArray([A00, A10, A10, A11], [2, 2])
};
})
eig = Eigen.ensure(invalidation, () => {
const eig = axes.covariance.selfAdjointEigenSolver(Eigen.ComputeEigenvectors);
const vals = eig.eigenvalues();
const vecs = eig.eigenvectors();
// Compute indices which sort the eigenvalues e.g. [1, 0] if they are reversed,
// and not to be confused with sorted eigenvalues themselves.
const eigOrder = argsort(Array.apply(null, vals.buffer));
// Return sorted eigenvalues and eigenvectors.
return {
vals: vals.map((i, j, val) => vals.get(eigOrder[i], j)),
vecs: vecs.map((i, j, val) => vecs.get(eigOrder[i], j))
};
})
viewof axisProjection = Eigen.ensure(invalidation, () =>
x
.sub(Matrixf.ones([n, 1]).mul(axes.centroid))
.mul(Matrixf.fromCol(eig.vecs.col(1)))
.show('(x - x_c) \\cdot v')
)
One platform to build and deploy the best data apps
Experiment and prototype by building visualizations in live JavaScript notebooks. Collaborate with your team and decide which concepts to build out.
Use Observable Framework to build data apps locally. Use data loaders to build in any language or library, including Python, SQL, and R.
Seamlessly deploy to Observable. Test before you ship, use automatic deploy-on-commit, and ensure your projects are always up-to-date.