Drawing cumulative density function (CDF).
- direct CDF (red line)
- via-histogram CDF (blue line)
Imspired by https://bl.ocks.org/jltran/bcd2a30fd9c08f8b9f87.
Last active
September 5, 2022 13:13
-
-
Save kkdd/99f4befb33eaff3d0aa31f99182cfd98 to your computer and use it in GitHub Desktop.
Drawing cumulative density function
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| <!DOCTYPE html> | |
| <!-- ref: https://bl.ocks.org/jltran/bcd2a30fd9c08f8b9f87 --> | |
| <meta charset="utf-8"> | |
| <style> | |
| svg { | |
| font: 10px sans-serif; | |
| } | |
| .axis path, .axis line { | |
| fill: none; | |
| stroke: #000; | |
| shape-rendering: crispEdges; | |
| } | |
| </style> | |
| <body> | |
| <script src="https://d3js.org/d3.v7.min.js"></script> | |
| <script> | |
| //Set dimensions | |
| const m = {top: 50, right: 50, bottom: 50, left: 50} | |
| , h = 500 - m.top - m.bottom | |
| , w = 960 - m.left - m.right | |
| , numBins = 50 | |
| , n = 10000; | |
| //Generate random numbers and calculate direct CDF | |
| const norm = d3.range(n).map(d3.randomNormal(0, 1)); | |
| norm.sort((first, second) => first - second); | |
| const cum = norm.slice(0, n - 1).map((w, i) => ({'x': w, 'y': (i+1)/n})); | |
| //Calculate histogram and derived CDR | |
| const xScale = d3.scaleLinear().domain([-3, 3]).range([0, w]); | |
| const data = d3.histogram().domain(xScale.domain()).thresholds(xScale.ticks(numBins))(norm).filter(d => d.x0 < d.x1); | |
| for(let i = 0, cum0 = 0; i < data.length; i++){ | |
| cum0 += data[i].length; | |
| data[i]['cum'] = cum0 / n; | |
| } | |
| //Scales and axes | |
| const yhist = d3.scaleLinear().range([h, 0]).domain([0, d3.max(data, d => d.length)]); | |
| const ycum = d3.scaleLinear().domain([0, 1]).range([h, 0]); | |
| const xAxis = d3.axisBottom(xScale); | |
| const yAxis = d3.axisLeft(yhist); | |
| const yAxis2 = d3.axisRight(ycum); | |
| //Draw svg | |
| const svg = d3.select("body").append("svg") | |
| .attr("width", w + m.left + m.right) | |
| .attr("height", h + m.top + m.bottom) | |
| .append("g") | |
| .attr("transform", "translate(" + m.left + "," + m.top + ")"); | |
| //Draw histogram | |
| svg.selectAll("rect") | |
| .data(data) | |
| .enter() | |
| .append("rect") | |
| .attr("x", 1) | |
| .attr("transform", d => "translate(" + xScale(d.x0) + "," + yhist(d.length) + ")") | |
| .attr("width", d => xScale(d.x1) - xScale(d.x0) - 1) | |
| .attr("height", d => h - yhist(d.length)) | |
| .style("fill", "#69b3a2") | |
| .style("opacity", 0.2); | |
| //Draw CDF's | |
| svg.append("path") | |
| .attr("d", d3.line().x(d => xScale(d.x1)).y(d => ycum(d.cum)).curve(d3.curveLinear)(data)) | |
| .attr("stroke", "#00ffff") | |
| .attr("stroke-width", 7) | |
| .attr("fill", "none"); | |
| svg.append("path") | |
| .attr("d", (d3.line().x(d => xScale(d.x)).y(d => ycum(erfinv2(d.y))).curve(d3.curveLinear)(cum))) | |
| .attr("stroke", "#ffc000") | |
| .attr("stroke-width", 1.2) | |
| .attr("fill", "none"); | |
| svg.append("path") | |
| .attr("d", (d3.line().x(d => xScale(d.x)).y(d => ycum(d.y)).curve(d3.curveLinear)(cum))) | |
| .attr("stroke", "#ff0000") | |
| .attr("stroke-width", 1.2) | |
| .attr("fill", "none"); | |
| //Draw axes | |
| svg.append("g") | |
| .attr("class", "x axis") | |
| .attr("transform", "translate(0," + h + ")") | |
| .call(xAxis); | |
| svg.append("g") | |
| .attr("class", "y axis") | |
| .call(yAxis) | |
| .append("text") | |
| .attr("transform", "rotate(-90)") | |
| .attr("y", 6) | |
| .attr("dy", ".71em") | |
| .style("text-anchor", "end") | |
| .text("Count (Histogram)"); | |
| svg.append("g") | |
| .attr("class", "y axis") | |
| .attr("transform", "translate(" + [w, 0] + ")") | |
| .call(yAxis2) | |
| .append("text") | |
| .attr("transform", "rotate(-90)") | |
| .attr("y", 4) | |
| .attr("dy", "-.71em") | |
| .style("text-anchor", "end") | |
| .text("CDF"); | |
| function erfinv2(x) { | |
| return (erfinvapprox(x * 2 - 1) * 2 / Math.sqrt(Math.PI) + 1 ) / 2; | |
| } | |
| function erfinvapprox(x) { | |
| // maximum relative error = .00013 | |
| const a = 0.147; | |
| const log1 = Math.log(1 - x**2); | |
| const b = 2 / (Math.PI * a) + log1/2 | |
| const sqrt1 = Math.sqrt(b**2 - log1/a) | |
| return Math.sqrt(sqrt1 - b) * Math.sign(x) | |
| } | |
| //if (x <= 0) { return -999 } else if (x >= 1) {return 999} | |
| </script> | |
| </body> |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment