weighted KDE & Voronoï maps

.block

license: gpl-3.0

README.md

This block experiments two things:

1. computation of Weighted KDE, an extension of standard KDE; in standard KDE, each datum counts for the same amount (i.e. 1), whereas in weighted KDE, each data counts for a specific amount (i.e. its weight); for example, one can make a standard KDE of the number of sales per day, or a weighted KDE of the total sales' profit per day
2. fill the weighted KDE curve/area with cells encoding data's weights (light weight -> small cell, heavy weight -> large cell); the objective is to give a sens of the underlying distribution that produces the weighted KDE; I use the d3-voronoï-map plugin to do so, but:
3. tweek it so that each site's x-coord remains unchanged durring the voronoï map computation (see file d"-voronoi-map-fixed-x.js)
4. define a specific initial positioning function (exactly encodes x-coord, and computes a random y-coord)

Usage : use the controller to hide/show objects, and hover a cell or a bin for details.

Indeed, the underlying dataset does not suit the experimentation. I have to find another one.

• forked from mbostock's block: Kernel Density Estimation
• see also http://bl.ocks.org/jfirebaugh/900762 and http://bl.ocks.org/z-m-k/5014368 for standard KDE computation/vizualisation with D3

==original README==

Kernel density estimation is a method of estimating the probability distribution of a random variable based on a random sample. In contrast to a histogram, kernel density estimation produces a smooth estimate. The smoothness can be tuned via the kernel’s bandwidth parameter. With the correct choice of bandwidth, important features of the distribution can be seen, while an incorrect choice results in undersmoothing or oversmoothing and obscured features.

This example shows a histogram and a kernel density estimation for times between eruptions of Old Faithful Geyser in Yellowstone National Park, taken from R’s faithful dataset. The data follow a bimodal distribution; short eruptions are followed by a wait time averaging about 55 minutes, and long eruptions by a wait time averaging about 80 minutes. In recent years, wait times have been increasing, possibly due to the effects of earthquakes on the geyser’s geohydrology.

This example is based on a Protovis version by John Firebaugh. See also a two-dimensional density estimation of this dataset using d3-contour.

d3-voronoi-map-fixed-x.js

// Code from github.com/Kcnarf/d3-voronoi-map/v1.2.0/build/d3-voronoi-map.js
// Adapted to fix x-coord of sites during the entire process
// look for 'adaptPosition' function and '//(fixed-x)' prefixed loc.



(function (global, factory) {
  typeof exports === 'object' && typeof module !== 'undefined' ? factory(exports, require('d3-polygon'), require('d3-weighted-voronoi')) :
  typeof define === 'function' && define.amd ? define(['exports', 'd3-polygon', 'd3-weighted-voronoi'], factory) :
  (factory((global.d3 = global.d3 || {}),global.d3,global.d3));
}(this, function (exports,d3Polygon,d3WeightedVoronoi) { 'use strict';

  function FlickeringMitigation () {
    /////// Inputs ///////
    this.growthChangesLength = DEFAULT_LENGTH;
    this.totalAvailableArea = NaN;
    
    //begin: internals
    this.lastAreaError = NaN;
    this.lastGrowth = NaN;
    this.growthChanges = [];
    this.growthChangeWeights = generateGrowthChangeWeights(this.growthChangesLength); //used to make recent changes weighter than older changes
    this.growthChangeWeightsSum = computeGrowthChangeWeightsSum(this.growthChangeWeights);
    //end: internals
  }

  var DEFAULT_LENGTH = 10;

  function direction(h0, h1) {
    return (h0 >= h1)? 1 : -1;
  }

  function generateGrowthChangeWeights(length) {
    var initialWeight = 3;   // a magic number
    var weightDecrement = 1; // a magic number
    var minWeight = 1;
    
    var weightedCount = initialWeight;
    var growthChangeWeights = [];
    
    for (var i=0; i<length; i++) {
      growthChangeWeights.push(weightedCount);
      weightedCount -= weightDecrement;
      if (weightedCount<minWeight) { weightedCount = minWeight; }
    }
    return growthChangeWeights;
  }

  function computeGrowthChangeWeightsSum (growthChangeWeights) {
    var growthChangeWeightsSum = 0;
    for (var i=0; i<growthChangeWeights.length; i++) {
      growthChangeWeightsSum += growthChangeWeights[i];
    }
    return growthChangeWeightsSum;
  }

  ///////////////////////
  ///////// API /////////
  ///////////////////////

  FlickeringMitigation.prototype.reset = function () {
    this.lastAreaError = NaN;
    this.lastGrowth = NaN;
    this.growthChanges = [];
    this.growthChangesLength = DEFAULT_LENGTH;
    this.growthChangeWeights = generateGrowthChangeWeights(this.growthChangesLength);
    this.growthChangeWeightsSum = computeGrowthChangeWeightsSum(this.growthChangeWeights);
    this.totalAvailableArea = NaN;
    
    return this;
  };

  FlickeringMitigation.prototype.clear = function () {
    this.lastAreaError = NaN;
    this.lastGrowth = NaN;
    this.growthChanges = [];
    
    return this;
  };

  FlickeringMitigation.prototype.length = function (_) {
    if (!arguments.length) { return this.growthChangesLength; }
    
    if (parseInt(_)>0) {
      this.growthChangesLength = Math.floor(parseInt(_));
      this.growthChangeWeights = generateGrowthChangeWeights(this.growthChangesLength);
      this.growthChangeWeightsSum = computeGrowthChangeWeightsSum(this.growthChangeWeights);
    } else {
      console.warn("FlickeringMitigation.length() accepts only positive integers; unable to handle "+_);
    }
    return this;
  };

  FlickeringMitigation.prototype.totalArea = function (_) {
    if (!arguments.length) { return this.totalAvailableArea; }
    
    if (parseFloat(_)>0) {
      this.totalAvailableArea = parseFloat(_);
    } else {
      console.warn("FlickeringMitigation.totalArea() accepts only positive numbers; unable to handle "+_);
    }
    return this;
  };

  FlickeringMitigation.prototype.add = function (areaError) {
    var secondToLastAreaError, secondToLastGrowth;

    secondToLastAreaError = this.lastAreaError;
    this.lastAreaError = areaError;
    if (!isNaN(secondToLastAreaError)) {
      secondToLastGrowth = this.lastGrowth;
      this.lastGrowth = direction(this.lastAreaError, secondToLastAreaError);
    }
    if (!isNaN(secondToLastGrowth)) {
      this.growthChanges.unshift(this.lastGrowth!=secondToLastGrowth);
    }

    if (this.growthChanges.length>this.growthChangesLength) {
      this.growthChanges.pop();
    }
    return this;
  };

  FlickeringMitigation.prototype.ratio = function () {
    var weightedChangeCount = 0;
    var ratio;

    if (this.growthChanges.length < this.growthChangesLength) { return 0; }
    if (this.lastAreaError > this.totalAvailableArea/10) { return 0; }

    for(var i=0; i<this.growthChangesLength; i++) {
      if (this.growthChanges[i]) {
        weightedChangeCount += this.growthChangeWeights[i];
      }
    }

    ratio = weightedChangeCount/this.growthChangeWeightsSum;

    if (ratio>0) {
      console.log("flickering mitigation ratio: "+Math.floor(ratio*1000)/1000);
    }

    return ratio;
  };

  function randomInitialPosition () {

    //begin: internals
    var clippingPolygon,
      extent,
      minX, maxX,
      minY, maxY,
      dx, dy;
    //end: internals

    ///////////////////////
    ///////// API /////////
    ///////////////////////

    function _random(d, i, arr, voronoiMap) {
      var shouldUpdateInternals = false;
      var x, y;

      if (clippingPolygon !== voronoiMap.clip()) {
        clippingPolygon = voronoiMap.clip();
        extent = voronoiMap.extent();
        shouldUpdateInternals = true;
      }

      if (shouldUpdateInternals) {
        updateInternals();
      }

      x = minX + dx * voronoiMap.prng()();
      y = minY + dy * voronoiMap.prng()();
      while (!d3Polygon.polygonContains(clippingPolygon, [x, y])) {
        x = minX + dx * voronoiMap.prng()();
        y = minY + dy * voronoiMap.prng()();
      }
      return [x, y];
    };

    ///////////////////////
    /////// Private ///////
    ///////////////////////

    function updateInternals() {
      minX = extent[0][0];
      maxX = extent[1][0];
      minY = extent[0][1];
      maxY = extent[1][1];
      dx = maxX - minX;
      dy = maxY - minY;
    };

    return _random;
  };

  function pie () {
    //begin: internals
    var startAngle = 0;
    var clippingPolygon,
      dataArray,
      dataArrayLength,
      clippingPolygonCentroid,
      halfIncircleRadius,
      angleBetweenData;
    //end: internals

    ///////////////////////
    ///////// API /////////
    ///////////////////////

    function _pie(d, i, arr, voronoiMap) {
      var shouldUpdateInternals = false;

      if (clippingPolygon !== voronoiMap.clip()) {
        clippingPolygon = voronoiMap.clip();
        shouldUpdateInternals |= true;
      }
      if (dataArray !== arr) {
        dataArray = arr;
        shouldUpdateInternals |= true;
      }

      if (shouldUpdateInternals) {
        updateInternals();
      }

      // add some randomness to prevent colinear/cocircular points
      // substract -0.5 so that the average jitter is still zero
      return [
        clippingPolygonCentroid[0] + Math.cos(startAngle + i * angleBetweenData) * halfIncircleRadius + (voronoiMap.prng()() - 0.5) * 1E-3,
        clippingPolygonCentroid[1] + Math.sin(startAngle + i * angleBetweenData) * halfIncircleRadius + (voronoiMap.prng()() - 0.5) * 1E-3
      ];
    };

    _pie.startAngle = function (_) {
      if (!arguments.length) {
        return startAngle;
      }

      startAngle = _;
      return _pie;
    };

    ///////////////////////
    /////// Private ///////
    ///////////////////////

    function updateInternals() {
      clippingPolygonCentroid = d3Polygon.polygonCentroid(clippingPolygon);
      halfIncircleRadius = computeMinDistFromEdges(clippingPolygonCentroid, clippingPolygon) / 2;
      dataArrayLength = dataArray.length;
      angleBetweenData = 2 * Math.PI / dataArrayLength;
    };

    function computeMinDistFromEdges(vertex, clippingPolygon) {
      var minDistFromEdges = Infinity,
        edgeIndex = 0,
        edgeVertex0 = clippingPolygon[clippingPolygon.length - 1],
        edgeVertex1 = clippingPolygon[edgeIndex];
      var distFromCurrentEdge;

      while (edgeIndex < clippingPolygon.length) {
        distFromCurrentEdge = vDistance(vertex, edgeVertex0, edgeVertex1);
        if (distFromCurrentEdge < minDistFromEdges) {
          minDistFromEdges = distFromCurrentEdge;
        }
        edgeIndex++;
        edgeVertex0 = edgeVertex1;
        edgeVertex1 = clippingPolygon[edgeIndex];
      }

      return minDistFromEdges;
    }

    //from https://stackoverflow.com/questions/849211/shortest-distance-between-a-point-and-a-line-segment
    function vDistance(vertex, edgeVertex0, edgeVertex1) {
      var x = vertex[0],
        y = vertex[1],
        x1 = edgeVertex0[0],
        y1 = edgeVertex0[1],
        x2 = edgeVertex1[0],
        y2 = edgeVertex1[1];
      var A = x - x1,
        B = y - y1,
        C = x2 - x1,
        D = y2 - y1;
      var dot = A * C + B * D;
      var len_sq = C * C + D * D;
      var param = -1;

      if (len_sq != 0) //in case of 0 length line
        param = dot / len_sq;

      var xx, yy;

      if (param < 0) { // this should not arise as clippingpolygon is convex
        xx = x1;
        yy = y1;
      } else if (param > 1) { // this should not arise as clippingpolygon is convex
        xx = x2;
        yy = y2;
      } else {
        xx = x1 + param * C;
        yy = y1 + param * D;
      }

      var dx = x - xx;
      var dy = y - yy;
      return Math.sqrt(dx * dx + dy * dy);
    }

    return _pie;
  }

  function halfAverageAreaInitialWeight () {
    //begin: internals
    var clippingPolygon,
      dataArray,
      siteCount,
      totalArea,
      halfAverageArea;
    //end: internals

    ///////////////////////
    ///////// API /////////
    ///////////////////////
    function _halfAverageArea(d, i, arr, voronoiMap) {
      var shouldUpdateInternals = false;
      if (clippingPolygon !== voronoiMap.clip()) {
        clippingPolygon = voronoiMap.clip();
        shouldUpdateInternals |= true;
      }
      if (dataArray !== arr) {
        dataArray = arr;
        shouldUpdateInternals |= true;
      }

      if (shouldUpdateInternals) {
        updateInternals();
      }

      return halfAverageArea;
    };

    ///////////////////////
    /////// Private ///////
    ///////////////////////

    function updateInternals() {
      siteCount = dataArray.length;
      totalArea = d3Polygon.polygonArea(clippingPolygon);
      halfAverageArea = totalArea / siteCount / 2; // half of the average area of the the clipping polygon
    }

    return _halfAverageArea;
  };

  function voronoiMap() {
    //begin: constants
    var DEFAULT_CONVERGENCE_RATIO = 0.01;
    var DEFAULT_MAX_ITERATION_COUNT = 50;
    var DEFAULT_MIN_WEIGHT_RATIO = 0.01;
    var DEFAULT_PRNG = Math.random;
    var DEFAULT_INITIAL_POSITION = randomInitialPosition();
    var DEFAULT_INITIAL_WEIGHT = halfAverageAreaInitialWeight();
    var RANDOM_INITIAL_POSITION = randomInitialPosition();
    var epsilon = 1;
    //end: constants

    /////// Inputs ///////
    var weight = function (d) {
      return d.weight;
    }; // accessor to the weight
    var convergenceRatio = DEFAULT_CONVERGENCE_RATIO; // targeted allowed error ratio; default 0.01 stops computation when cell areas error <= 1% clipping polygon's area
    var maxIterationCount = DEFAULT_MAX_ITERATION_COUNT; // maximum allowed iteration; stops computation even if convergence is not reached; use a large amount for a sole converge-based computation stop
    var minWeightRatio = DEFAULT_MIN_WEIGHT_RATIO; // used to compute the minimum allowed weight; default 0.01 means 1% of max weight; handle near-zero weights, and leaves enought space for cell hovering
    var prng = DEFAULT_PRNG; // pseudorandom number generator
    var initialPosition = DEFAULT_INITIAL_POSITION; // accessor to the initial position; defaults to a random position inside the clipping polygon
    var initialWeight = DEFAULT_INITIAL_WEIGHT; // accessor to the initial weight; defaults to the average area of the clipping polygon
    var tick = function (polygons, i) {
      return true;
    }; // hook called at each iteration's end (i = iteration count)

    //begin: internals
    var weightedVoronoi = d3WeightedVoronoi.weightedVoronoi();
    var siteCount,
      totalArea,
      areaErrorTreshold,
      flickeringMitigation = new FlickeringMitigation();
    //end: internals

    //begin: algorithm conf.
    var handleOverweightedVariant = 1; // this option still exists 'cause for further experiments
    var handleOverweighted;
    //end: algorithm conf.

    //begin: utils
    function sqr(d) {
      return Math.pow(d, 2);
    }

    function squaredDistance(s0, s1) {
      return sqr(s1.x - s0.x) + sqr(s1.y - s0.y);
    }
    //end: utils

    ///////////////////////
    ///////// API /////////
    ///////////////////////

    function _voronoiMap(data) {
      //begin: handle algorithm's variants
      setHandleOverweighted();
      //end: handle algorithm's variants

      siteCount = data.length;
      (totalArea = Math.abs(d3Polygon.polygonArea(weightedVoronoi.clip()))), (areaErrorTreshold = convergenceRatio * totalArea);
      flickeringMitigation.clear().totalArea(totalArea);

      var iterationCount = 0,
        polygons = initialize(data),
        converged = false;
      var areaError;

      tick(polygons, iterationCount);

      while (!(converged || iterationCount >= maxIterationCount)) {
        polygons = adapt(polygons, flickeringMitigation.ratio());
        iterationCount++;
        areaError = computeAreaError(polygons);
        flickeringMitigation.add(areaError);
        converged = areaError < areaErrorTreshold;
        // console.log("error %: "+Math.round(areaError*100*1000/totalArea)/1000);
        tick(polygons, iterationCount);
      }

      return {
        polygons: polygons,
        iterationCount: iterationCount,
        convergenceRatio: areaError / totalArea
      };
    }

    _voronoiMap.weight = function (_) {
      if (!arguments.length) {
        return weight;
      }

      weight = _;
      return _voronoiMap;
    };

    _voronoiMap.convergenceRatio = function (_) {
      if (!arguments.length) {
        return convergenceRatio;
      }

      convergenceRatio = _;
      return _voronoiMap;
    };

    _voronoiMap.maxIterationCount = function (_) {
      if (!arguments.length) {
        return maxIterationCount;
      }

      maxIterationCount = _;
      return _voronoiMap;
    };

    _voronoiMap.minWeightRatio = function (_) {
      if (!arguments.length) {
        return minWeightRatio;
      }

      minWeightRatio = _;
      return _voronoiMap;
    };

    _voronoiMap.tick = function (_) {
      if (!arguments.length) {
        return tick;
      }

      tick = _;
      return _voronoiMap;
    };

    _voronoiMap.clip = function (_) {
      if (!arguments.length) {
        return weightedVoronoi.clip();
      }

      weightedVoronoi.clip(_);
      return _voronoiMap;
    };

    _voronoiMap.extent = function (_) {
      if (!arguments.length) {
        return weightedVoronoi.extent();
      }

      weightedVoronoi.extent(_);
      return _voronoiMap;
    };

    _voronoiMap.size = function (_) {
      if (!arguments.length) {
        return weightedVoronoi.size();
      }

      weightedVoronoi.size(_);
      return _voronoiMap;
    };

    _voronoiMap.prng = function (_) {
      if (!arguments.length) {
        return prng;
      }

      prng = _;
      return _voronoiMap;
    };

    _voronoiMap.initialPosition = function (_) {
      if (!arguments.length) {
        return initialPosition;
      }

      initialPosition = _;
      return _voronoiMap;
    };

    _voronoiMap.initialWeight = function (_) {
      if (!arguments.length) {
        return initialWeight;
      }

      initialWeight = _;
      return _voronoiMap;
    };

    ///////////////////////
    /////// Private ///////
    ///////////////////////

    function adapt(polygons, flickeringMitigationRatio) {
      var adaptedMapPoints;

      adaptPositions(polygons, flickeringMitigationRatio);
      adaptedMapPoints = polygons.map(function (p) {
        return p.site.originalObject;
      });
      polygons = weightedVoronoi(adaptedMapPoints);
      if (polygons.length < siteCount) {
        console.log('at least 1 site has no area, which is not supposed to arise');
        debugger;
      }

      adaptWeights(polygons, flickeringMitigationRatio);
      adaptedMapPoints = polygons.map(function (p) {
        return p.site.originalObject;
      });
      polygons = weightedVoronoi(adaptedMapPoints);
      if (polygons.length < siteCount) {
        console.log('at least 1 site has no area, which is not supposed to arise');
        debugger;
      }

      return polygons;
    }

    function adaptPositions(polygons, flickeringMitigationRatio) {
      var newMapPoints = [],
        flickeringInfluence = 0.5;
      var flickeringMitigation, d, polygon, mapPoint, centroid, dx, dy;

      flickeringMitigation = flickeringInfluence * flickeringMitigationRatio;
      d = 1 - flickeringMitigation; // in [0.5, 1]
      for (var i = 0; i < siteCount; i++) {
        polygon = polygons[i];
        mapPoint = polygon.site.originalObject;
        centroid = d3Polygon.polygonCentroid(polygon);

        //(fixed-x)dx = centroid[0] - mapPoint.x;
        dy = centroid[1] - mapPoint.y;

        //begin: handle excessive change;
        //(fixed-x)dx *= d;
        dy *= d;
        //end: handle excessive change;

        //(fixed-x)mapPoint.x += dx;
        mapPoint.y += dy;

        newMapPoints.push(mapPoint);
      }

      handleOverweighted(newMapPoints);
    }

    function adaptWeights(polygons, flickeringMitigationRatio) {
      var newMapPoints = [],
        flickeringInfluence = 0.1;
      var flickeringMitigation, polygon, mapPoint, currentArea, adaptRatio, adaptedWeight;

      flickeringMitigation = flickeringInfluence * flickeringMitigationRatio;
      for (var i = 0; i < siteCount; i++) {
        polygon = polygons[i];
        mapPoint = polygon.site.originalObject;
        currentArea = d3Polygon.polygonArea(polygon);
        adaptRatio = mapPoint.targetedArea / currentArea;

        //begin: handle excessive change;
        adaptRatio = Math.max(adaptRatio, 1 - flickeringInfluence + flickeringMitigation); // in [(1-flickeringInfluence), 1]
        adaptRatio = Math.min(adaptRatio, 1 + flickeringInfluence - flickeringMitigation); // in [1, (1+flickeringInfluence)]
        //end: handle excessive change;

        adaptedWeight = mapPoint.weight * adaptRatio;
        adaptedWeight = Math.max(adaptedWeight, epsilon);

        mapPoint.weight = adaptedWeight;

        newMapPoints.push(mapPoint);
      }

      handleOverweighted(newMapPoints);
    }

    // heuristics: lower heavy weights
    function handleOverweighted0(mapPoints) {
      var fixCount = 0;
      var fixApplied, tpi, tpj, weightest, lightest, sqrD, adaptedWeight;
      do {
        fixApplied = false;
        for (var i = 0; i < siteCount; i++) {
          tpi = mapPoints[i];
          for (var j = i + 1; j < siteCount; j++) {
            tpj = mapPoints[j];
            if (tpi.weight > tpj.weight) {
              weightest = tpi;
              lightest = tpj;
            } else {
              weightest = tpj;
              lightest = tpi;
            }
            sqrD = squaredDistance(tpi, tpj);
            if (sqrD < weightest.weight - lightest.weight) {
              // adaptedWeight = sqrD - epsilon; // as in ArlindNocaj/Voronoi-Treemap-Library
              // adaptedWeight = sqrD + lightest.weight - epsilon; // works, but below heuristics performs better (less flickering)
              adaptedWeight = sqrD + lightest.weight / 2;
              adaptedWeight = Math.max(adaptedWeight, epsilon);
              weightest.weight = adaptedWeight;
              fixApplied = true;
              fixCount++;
              break;
            }
          }
          if (fixApplied) {
            break;
          }
        }
      } while (fixApplied);

      /*
      if (fixCount>0) {
        console.log("# fix: "+fixCount);
      }
      */
    }

    // heuristics: increase light weights
    function handleOverweighted1(mapPoints) {
      var fixCount = 0;
      var fixApplied, tpi, tpj, weightest, lightest, sqrD, overweight;
      do {
        fixApplied = false;
        for (var i = 0; i < siteCount; i++) {
          tpi = mapPoints[i];
          for (var j = i + 1; j < siteCount; j++) {
            tpj = mapPoints[j];
            if (tpi.weight > tpj.weight) {
              weightest = tpi;
              lightest = tpj;
            } else {
              weightest = tpj;
              lightest = tpi;
            }
            sqrD = squaredDistance(tpi, tpj);
            if (sqrD < weightest.weight - lightest.weight) {
              overweight = weightest.weight - lightest.weight - sqrD;
              lightest.weight += overweight + epsilon;
              fixApplied = true;
              fixCount++;
              break;
            }
          }
          if (fixApplied) {
            break;
          }
        }
      } while (fixApplied);

      /*
      if (fixCount>0) {
        console.log("# fix: "+fixCount);
      }
      */
    }

    function computeAreaError(polygons) {
      //convergence based on summation of all sites current areas
      var areaErrorSum = 0;
      var polygon, mapPoint, currentArea;
      for (var i = 0; i < siteCount; i++) {
        polygon = polygons[i];
        mapPoint = polygon.site.originalObject;
        currentArea = d3Polygon.polygonArea(polygon);
        areaErrorSum += Math.abs(mapPoint.targetedArea - currentArea);
      }
      return areaErrorSum;
    }

    function setHandleOverweighted() {
      switch (handleOverweightedVariant) {
        case 0:
          handleOverweighted = handleOverweighted0;
          break;
        case 1:
          handleOverweighted = handleOverweighted1;
          break;
        default:
          console.log("Variant of 'handleOverweighted' is unknown");
      }
    }

    function initialize(data) {
      var maxWeight = data.reduce(function (max, d) {
          return Math.max(max, weight(d));
        }, -Infinity),
        minAllowedWeight = maxWeight * minWeightRatio;
      var weights, mapPoints;

      //begin: extract weights
      weights = data.map(function (d, i, arr) {
        return {
          index: i,
          weight: Math.max(weight(d), minAllowedWeight),
          initialPosition: initialPosition(d, i, arr, _voronoiMap),
          initialWeight: initialWeight(d, i, arr, _voronoiMap),
          originalData: d
        };
      });
      //end: extract weights

      // create map-related points
      // (with targetedArea, initial position and initialWeight)
      mapPoints = createMapPoints(weights);
      return weightedVoronoi(mapPoints);
    }

    function createMapPoints(basePoints) {
      var totalWeight = basePoints.reduce(function (acc, bp) {
        return (acc += bp.weight);
      }, 0);
      var initialPosition;

      return basePoints.map(function (bp, i, bps) {
        initialPosition = bp.initialPosition;

        if (!d3Polygon.polygonContains(weightedVoronoi.clip(), initialPosition)) {
          initialPosition = RANDOM_INITIAL_POSITION(bp, i, bps, _voronoiMap);
        }

        return {
          index: bp.index,
          targetedArea: totalArea * bp.weight / totalWeight,
          data: bp,
          x: initialPosition[0],
          y: initialPosition[1],
          weight: bp.initialWeight // ArlindNocaj/Voronoi-Treemap-Library uses an epsilonesque initial weight; using heavier initial weights allows faster weight adjustements, hence faster stabilization
        };
      });
    }

    return _voronoiMap;
  }

  exports.voronoiMap = voronoiMap;
  exports.voronoiMapInitialPositionRandom = randomInitialPosition;
  exports.voronoiMapInitialPositionPie = pie;

  Object.defineProperty(exports, '__esModule', { value: true });

}));

faithful.json

[79,54,74,62,85,55,88,85,51,85,54,84,78,47,83,52,62,84,52,79,51,47,78,69,74,83,55,76,78,79,73,77,66,80,74,52,48,80,59,90,80,58,84,58,73,83,64,53,82,59,75,90,54,80,54,83,71,64,77,81,59,84,48,82,60,92,78,78,65,73,82,56,79,71,62,76,60,78,76,83,75,82,70,65,73,88,76,80,48,86,60,90,50,78,63,72,84,75,51,82,62,88,49,83,81,47,84,52,86,81,75,59,89,79,59,81,50,85,59,87,53,69,77,56,88,81,45,82,55,90,45,83,56,89,46,82,51,86,53,79,81,60,82,77,76,59,80,49,96,53,77,77,65,81,71,70,81,93,53,89,45,86,58,78,66,76,63,88,52,93,49,57,77,68,81,81,73,50,85,74,55,77,83,83,51,78,84,46,83,55,81,57,76,84,77,81,87,77,51,78,60,82,91,53,78,46,77,84,49,83,71,80,49,75,64,76,53,94,55,76,50,82,54,75,78,79,78,78,70,79,70,54,86,50,90,54,54,77,79,64,75,47,86,63,85,82,57,82,67,74,54,83,73,73,88,80,71,83,56,79,78,84,58,83,43,60,75,81,46,90,46,74]

index.html

<!DOCTYPE html>
<html>
  <head>
    <meta charset="utf-8">
    <meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
    <title>weighted KDE + Voronoï maps</title>
    <meta name="description" content="Weighted KDE filled with Voronoï maps in D3.js">
    <style>
      .hide {
        display: none;
      }

      .standard {
        fill: #bbb;
        stroke: #bbb;
      }

      .weighted {
        fill: #ddd;
        stroke: #ddd;
      }

      .bin {
        stroke: white;
      }

      .kde {
        stroke-width: 1.5;
        stroke-linejoin: round;
        fill: none;
        stroke: black;
      }

      text.standard {
        stroke: none;
      }
      text.weighted {
        stroke: none;
      }

      .cell {
        fill: transparent;
        stroke-width: 1;
        stroke: grey;
      }

      .axis--y .domain {
        display: none;
      }
      
      .dg .property-name {
        width: 80% !important;
      }
      .dg .c {
        width: 20% !important;
      }
      
    </style>
  </head>
  <body>
    <svg width="960" height="500"></svg>
    <script src="https://d3js.org/d3.v5.min.js"></script>
    <script src="https://rawcdn.githack.com/Kcnarf/d3-weighted-voronoi/v1.0.0/build/d3-weighted-voronoi.js"></script>
    <script src="d3-voronoi-map-fixed-x.js"></script>
    <script src="//cdnjs.cloudflare.com/ajax/libs/seedrandom/2.4.3/seedrandom.min.js"></script>
    <script src="//cdnjs.cloudflare.com/ajax/libs/dat-gui/0.6.5/dat.gui.min.js"></script>
    <script>

      var settings = {
        standard: {
          showBins: true,
          showKDE: true
        },
        weighted: {
        	showWeightedBins: true,
        	showWeightedKDE: true,
        	showVoronoiMaps: true
      	}
      };
      var binCount = 40;

      var svg = d3.select("svg"),
          width = +svg.attr("width"),
          height = +svg.attr("height"),
          margin = {top: 20, right: 30, bottom: 30, left: 40};

      var x = d3.scaleLinear()
          .domain([30, 110])
          .range([margin.left, width - margin.right]);

      var y = d3.scaleLinear()
          .range([height - margin.bottom, margin.top]);

      var valueAccessor = d => d.value,
          weightAccessor = d => d.weight;

      var bins,							// will store histograms data
          density,					// will store KDE data
          weightedDensity,	// will store weighted KDE data
          maps;							// will store each Voronoi map's data of each bin

      d3.json("faithful.json").then(function(data) {
        // I know the data is not suitable for such an experimentation :-(
        // Data dealing with sales (with the possibility to count the number the sales, or to count the total sales' prices), would be much more better
        faithful = data.map((d)=>{
          return {
            value: d,
            weight: Math.round(d/10)	// weird weight feature :-(
          }
        });

        //begin: compute histograms data
        bins = d3.histogram()
                  .domain(x.domain())
                  .value(valueAccessor)
                  .thresholds(binCount)
              (faithful);
        bins.forEach(b=>{
          b.totalWeight = d3.sum(b, weightAccessor);
        })
        //end: compute histograms data

        y.domain([0, d3.max(bins, b=>b.totalWeight)*1.1]);

        //compute KDE's data, using a weight of 1 for each datum
        density = kernelDensityEstimator(x.ticks(binCount), scaledEpanechnikov(7, valueAccessor), d=>1)(faithful);

        //compute weighted KDE's data, using a particular weight for each datum
        weightedDensity = kernelDensityEstimator(x.ticks(binCount), scaledEpanechnikov(7, valueAccessor), weightAccessor)(faithful);

        //begin: weighted voronoï map for each right trapeze bin
        var seededprng = new Math.seedrandom('my seed'),
            voronoiMap = d3.voronoiMap().weight(weightAccessor).prng(seededprng);
        maps = [];
        weightedDensity.forEach(function (d,i) {
          if (i > weightedDensity.length-2) return;

          //begin: retrieve weighted KDE of bin's extent for right trapeze bins
          var e = weightedDensity[i+1],
              x0 = x(d[0]),
              y0 = y(d[1]),
              x1 = x(e[0]),
              y1 = y(e[1]),
              baseY = y(0),
              data = bins[i];
          //begin: retrieve weighted KDE of bin's extent for right trapeze bins
          var map;

          if (data.length <= 0) return;	// no real data

          //begin: voronoï map computation of a specific bin, considering weighted KDE
          map = voronoiMap
                  .clip([[x0,baseY], [x0,y0], [x1,y1], [x1,baseY]])
                  .weight(weightAccessor)
          				.initialPosition(randomYInitialPosition())
                (data).polygons;
          //begin: voronoï map computation of a specific bin, considering weighted KDE

          maps.push(map);
        })
        //end: weighted voronoï map for each right trapeze bin

        drawBins();
        drawVoronoiMaps();
        drawKDEs();
        drawAxes();
        drawDatGui();
      });


      /******************************/
      /*        Weighted KDE        */
      /******************************/

      function kernelDensityEstimator(samples, kernel, weightAccessor) {
        return function(data) {
          return samples.map(function(sample) {
            return [sample, d3.sum(data, function(d) {
              //return a value depending on kernel function AND weight
              return kernel(sample, d) * weightAccessor(d);
            })];
          });
        };
      }

      function scaledEpanechnikov(bw, valueAccessor) {
        // 'bd' states for 'bandwidth'
        return function(sample, d) {
          var diff = sample - valueAccessor(d);
          return epanechnikov(diff/bw)/bw;
        };
      }

      function epanechnikov(diff) {
        return Math.abs(diff) <= 1 ? 0.75 * (1 - diff * diff) : 0;
      };
      
      /******************************************************/
      /* Voronoï map initial positioning with fixed x-coord */
      /******************************************************/
      
      function randomYInitialPosition () {
        var clippingPolygon,
          extent,
          minY, maxY,
          dy;
        
        function _randomY(d, i, arr, voronoiMap) {
          var shouldUpdateInternals = false;
          var y;

          if (clippingPolygon !== voronoiMap.clip()) {
            clippingPolygon = voronoiMap.clip();
            extent = voronoiMap.extent();
            shouldUpdateInternals = true;
          }

          if (shouldUpdateInternals) {
            updateInternals();
          }

          y = minY + dy * voronoiMap.prng()();
          while (!d3.polygonContains(clippingPolygon, [x(valueAccessor(d)), y])) {
            y = minY + dy * voronoiMap.prng()();
          }
          return [x(valueAccessor(d)), y];
        };

        ///////////////////////
        /////// Private ///////
        ///////////////////////

        function updateInternals() {
          minY = extent[0][1];
          maxY = extent[1][1];
          dy = maxY - minY;
        };

        return _randomY;
      };


      /******************************/
      /*          Drawings          */
      /******************************/

      function drawBins() {
        //begin: weighted bin
        svg.append("g")
        		.attr("id", "weighted-histogram")
            .classed("hide", !settings.weighted.showWeightedBins)
            .selectAll("rect")
            .data(bins)
            .enter().append("rect")
              .classed("weighted bin", true)
              .attr("x", function(d) { return x(d.x0); })
              .attr("y", function(d) { return y(d.totalWeight); })
              .attr("width", function(d) { return x(d.x1) - x(d.x0); })
              .attr("height", function(d) { return y(0) - y(d.totalWeight); })
        			.append("title")
            		.text(d => { return "weighted count: "+d.totalWeight+"\ncount: "+d.length; });;
        //end: weighted bin

        //begin: standard bin
        svg.append("g")
        		.attr("id", "standard-histogram")
            .classed("hide", !settings.standard.showBins)
            .selectAll("rect")
            .data(bins)
            .enter().append("rect")
              .classed("standard bin", true)
              .attr("x", function(d) { return x(d.x0); })
              .attr("y", function(d) { return y(d.length); })
              .attr("width", function(d) { return x(d.x1) - x(d.x0); })
              .attr("height", function(d) { return y(0) - y(d.length); })
        			.append("title")
            		.text(d => { return "count: " + d.length; });
        //end: standard bin
      }

      function drawKDEs () {
        //begin: weighted KDE
        svg.append("path")
        		.attr("id", "weighted-kde")
            .classed("weighted kde", true)
            .classed("hide", !settings.weighted.showWeightedKDE)
            .datum(weightedDensity)
            .attr("d",  d3.line()
                .curve(d3.curveBasis)
                .x(function(d) { return x(d[0]); })
                .y(function(d) { return y(d[1]); }));
        //begin: weighted KDE

        //begin: standard KDE
        svg.append("path")
        		.attr("id", "standard-kde")
            .classed("standard kde", true)
            .classed("hide", !settings.standard.showKDE)
            .datum(density)
            .attr("d",  d3.line()
                .curve(d3.curveLinear)
                .x(function(d) { return x(d[0]); })
                .y(function(d) { return y(d[1]); }));
        //end: standard KDE
      }

      function drawVoronoiMaps () {
        //begin: draw each cells of each voronoï map
        var mapContainer = svg.append("g").attr("id", "voronoi-map-container").classed("hide", !settings.weighted.showVoronoiMaps);
        maps.forEach(function(d){
          mapContainer.append("g").classed("map", true)
            .selectAll("path")
            .data(d)
            .enter()
              .append("path")
              .classed("cell", true)
              .attr("d", d => { return d3.line().curve(d3.curveLinear)(d) + "z"; })
          		.append("title")
            		.text(d => { return "weight: " + weightAccessor(d.site.originalObject.data.originalData); })
        })
        //end: draw each cells of each voronoï map
      }

      function drawAxes () {
        svg.append("g")
          .attr("class", "axis axis--x")
          .attr("transform", "translate(0," + (height - margin.bottom) + ")")
          .call(d3.axisBottom(x))
          .append("text")
            .attr("x", width - margin.right)
            .attr("y", -6)
            .attr("fill", "#000")
            .attr("text-anchor", "end")
            .attr("font-weight", "bold")
            .text("Time between eruptions (min.)");

        var yAxis = svg.append("g")
          .attr("class", "axis axis--y")
          .attr("transform", "translate(" + margin.left + ",0)")
          .call(d3.axisLeft(y));
        yAxis.append("text")
            .classed("standard", true)
            .attr("x", -margin.top-80)
            .attr("y", 10)
            .style("text-anchor", "end")
            .style("font-weight", "bold")
            .style("transform", "rotate(-90deg)")
            .text("Count (i.e. number of eruptions) /");
        yAxis.append("text")
            .classed("weighted", true)
            .attr("x", -margin.top)
            .attr("y", 10)
            .style("text-anchor", "end")
            .style("font-weight", "bold")
            .style("transform", "rotate(-90deg)")
            .text("Weighted Count");
      }
      
      function drawDatGui() {
        var controls = new dat.GUI({width: 200});
        var standardControler = controls.addFolder("Standard");
        standardControler.add(settings.standard, "showBins")
          .listen().onChange(function(value) {
            d3.select("#standard-histogram").classed("hide", !value);
          });
        standardControler.add(settings.standard, "showKDE")
          .listen().onChange(function(value) {
            d3.select("#standard-kde").classed("hide", !value);
          });
        standardControler.open();
        var weightedControler = controls.addFolder("Weighted");
        weightedControler.add(settings.weighted, "showWeightedBins")
          .listen().onChange(function(value) {
            d3.select("#weighted-histogram").classed("hide", !value);
          });
        weightedControler.add(settings.weighted, "showWeightedKDE")
          .listen().onChange(function(value) {
            d3.select("#weighted-kde").classed("hide", !value);
          });
        weightedControler.add(settings.weighted, "showVoronoiMaps")
          .listen().onChange(function(value) {
            d3.select("#voronoi-map-container").classed("hide", !value);
          });
        weightedControler.open();
      }

    </script>
  </body>
</html>

thumbnail.png