31:T4e8,{ var i = polynome_degree; return tex` f(x) = ${polynome_terms[0].toFixed(8)} ${polynome_terms[1] >0 ? '+':''} ${polynome_terms[1].toFixed(10)} \, x ${polynome_terms[2] >0 && i >=2 ? '+':''} ${i >=2 ? `${polynome_terms[2].toFixed(10)} \\, x^2`:''} ${polynome_terms[3] >0 && i >=3 ? '+':''} ${i >=3 ? `${polynome_terms[3].toFixed(10)} \\, x^3`:''} ${polynome_terms[4] >0 && i >=4 ? '+':''} ${i >=4 ? `${polynome_terms[4].toFixed(10)} \\, x^4`:''} ${polynome_terms[5] >0 && i >=5 ? '+':''} ${i >=5 ? `${polynome_terms[5].toFixed(12)} \\, x^5`:''} ${i >=6 ? ' \\newline \\qquad \\quad \\;':''} ${polynome_terms[6] >0 && i >=6 ? '+':''} ${i >=6 ? `${polynome_terms[6].toFixed(15)} \\, x^6`:''} ${polynome_terms[7] >0 && i >=7 ? '+':''} ${i >=7 ? `${polynome_terms[7].toFixed(18)} \\, x^7`:''} ${polynome_terms[8] >0 && i >=8 ? '+':''} ${i >=8 ? `${polynome_terms[8].toFixed(18)} \\, x^8`:''} ${i >=9 ? ' \\newline \\qquad \\quad \\;':''} ${polynome_terms[9] >0 && i >=9 ? '+':''} ${i >=9 ? `${polynome_terms[9].toFixed(20)} \\, x^9`:''} ${polynome_terms[10] >0 && i >=10 ? '+':''} ${i >=10 ? `${polynome_terms[10].toFixed(24)} \\, x^{10}`:''} ` }32:T224e,// Source: https://github.com/rfink/polyfit.js // See also: https://arachnoid.com/polysolve/ Polyfit = (function () { /** * Polyfit * @constructor * @param {number[]|Float32Array|Float64Array} x * @param {number[]|Float32Array|Float64Array} y */ function Polyfit(x, y) { // Make sure we return an instance if (!(this instanceof Polyfit)) { return new Polyfit(x, y); } // Check that x any y are both arrays of the same type if (!((x instanceof Array && y instanceof Array) || (x instanceof Float32Array && y instanceof Float32Array) || (x instanceof Float64Array && y instanceof Float64Array))) { throw new Error('x and y must be arrays'); } if (x instanceof Float32Array) { this.FloatXArray = Float32Array; } else if (x instanceof Float64Array) { this.FloatXArray = Float64Array; } // Make sure we have equal lengths if (x.length !== y.length) { throw new Error('x and y must have the same length'); } this.x = x; this.y = y; } /** * Perform gauss-jordan division * * @param {number[][]|Float32Array[]|Float64Array[]} matrix - gets modified * @param {number} row * @param {number} col * @param {number} numCols * @returns void */ Polyfit.gaussJordanDivide = function (matrix, row, col, numCols) { for (var i = col + 1; i < numCols; i++) { matrix[row][i] /= matrix[row][col]; } matrix[row][col] = 1; }; /** * Perform gauss-jordan elimination * * @param {number[][]|Float64Array[]} matrix - gets modified * @param {number} row * @param {number} col * @param {number} numRows * @param {number} numCols * @returns void */ Polyfit.gaussJordanEliminate = function (matrix, row, col, numRows, numCols) { for (var i = 0; i < numRows; i++) { if (i !== row && matrix[i][col] !== 0) { for (var j = col + 1; j < numCols; j++) { matrix[i][j] -= matrix[i][col] * matrix[row][j]; } matrix[i][col] = 0; } } }; /** * Perform gauss-jordan echelon method * * @param {number[][]|Float32Array[]|Float64Array[]} matrix - gets modified * @returns {number[][]|Float32Array[]|Float64Array[]} matrix */ Polyfit.gaussJordanEchelonize = function (matrix) { var rows = matrix.length; var cols = matrix[0].length; var i = 0; var j = 0; var k; var swap; while (i < rows && j < cols) { k = i; // Look for non-zero entries in col j at or below row i while (k < rows && matrix[k][j] === 0) { k++; } // If an entry is found at row k if (k < rows) { // If k is not i, then swap row i with row k if (k !== i) { swap = matrix[i]; matrix[i] = matrix[k]; matrix[k] = swap; } // If matrix[i][j] is != 1, divide row i by matrix[i][j] if (matrix[i][j] !== 1) { Polyfit.gaussJordanDivide(matrix, i, j, cols); } // Eliminate all other non-zero entries Polyfit.gaussJordanEliminate(matrix, i, j, rows, cols); i++; } j++; } return matrix; }; /** * Perform regression * * @param {number} x * @param {number[]|Float32Array[]|Float64Array[]} terms * @returns {number} */ Polyfit.regress = function (x, terms) { var a = 0; var exp = 0; for (var i = 0, len = terms.length; i < len; i++) { a += terms[i] * Math.pow(x, exp++); } return a; }; /** * Compute correlation coefficient * * @param {number[]|Float32Array[]|Float64Array[]} terms * @returns {number} */ Polyfit.prototype.correlationCoefficient = function (terms) { var r = 0; var n = this.x.length; var sx = 0; var sx2 = 0; var sy = 0; var sy2 = 0; var sxy = 0; var x; var y; for (var i = 0; i < n; i++) { x = Polyfit.regress(this.x[i], terms); y = this.y[i]; sx += x; sy += y; sxy += x * y; sx2 += x * x; sy2 += y * y; } var div = Math.sqrt((sx2 - (sx * sx) / n) * (sy2 - (sy * sy) / n)); if (div !== 0) { r = Math.pow((sxy - (sx * sy) / n) / div, 2); } return r; }; /** * Run standard error function * * @param {number[]|Float32Array[]|Float64Array[]} terms * @returns number */ Polyfit.prototype.standardError = function (terms) { var r = 0; var n = this.x.length; if (n > 2) { var a = 0; for (var i = 0; i < n; i++) { a += Math.pow((Polyfit.regress(this.x[i], terms) - this.y[i]), 2); } r = Math.sqrt(a / (n - 2)); } return r; }; /** * Compute coefficients for given data matrix * * @param {number} p * @returns {number[]|Float32Array|Float64Array} */ Polyfit.prototype.computeCoefficients = function (p) { var n = this.x.length; var r; var c; var rs = 2 * (++p) - 1; var i; var m = []; // Initialize array with 0 values if (this.FloatXArray) { // fast FloatXArray-Matrix init var bytesPerRow = (p + 1) * this.FloatXArray.BYTES_PER_ELEMENT; var buffer = new ArrayBuffer(p * bytesPerRow); for (i = 0; i < p; i++) { m[i] = new this.FloatXArray(buffer, i * bytesPerRow, p + 1); } } else { var zeroRow = []; for (i = 0; i <= p; i++) { zeroRow[i] = 0; } m[0] = zeroRow; for (i = 1; i < p; i++) { // copy zeroRow m[i] = zeroRow.slice(); } } var mpc = [n]; for (i = 1; i < rs; i++) { mpc[i] = 0; } for (i = 0; i < n; i++) { var x = this.x[i]; var y = this.y[i]; // Process precalculation array for (r = 1; r < rs; r++) { mpc[r] += Math.pow(x, r); } // Process RH column cells m[0][p] += y; for (r = 1; r < p; r++) { m[r][p] += Math.pow(x, r) * y; } } // Populate square matrix section for (r = 0; r < p; r++) { for (c = 0; c < p; c++) { m[r][c] = mpc[r + c]; } } Polyfit.gaussJordanEchelonize(m); var terms = this.FloatXArray && new this.FloatXArray(m.length) || []; for (i = m.length - 1; i >= 0; i--) { terms[i] = m[i][p]; } return terms; }; /** * Using given degree of fitment, return a function that will calculate * the y for a given x * * @param {number} degree > 0 * @returns {Function} f(x) = */ Polyfit.prototype.getPolynomial = function (degree) { if (isNaN(degree) || degree < 0) { throw new Error('Degree must be a positive integer'); } var terms = this.computeCoefficients(degree); var eqParts = []; eqParts.push(terms[0].toPrecision()); for (var i = 1, len = terms.length; i < len; i++) { eqParts.push(terms[i] + ' * Math.pow(x, ' + i + ')'); } var expr = 'return ' + eqParts.join(' + ') + ';'; /* jshint evil: true */ return new Function('x', expr); /* jshint evil: false */ }; /** * Convert the polynomial to a string expression, mostly useful for visual * debugging * * @param {number} degree * @returns {string} */ Polyfit.prototype.toExpression = function (degree) { if (isNaN(degree) || degree < 0) { throw new Error('Degree must be a positive integer'); } var terms = this.computeCoefficients(degree); var eqParts = []; var len = terms.length; eqParts.push(terms[0].toPrecision()); for (var i = 1; i < len; i++) { eqParts.push(terms[i] + 'x^' + i); } return eqParts.join(' + '); }; return Polyfit; })();33:[] 2f:{"title":"Polynomial regression in JavaScript","cells":[{"id":0,"value":"# Polynomial regression in JavaScript","pinned":false,"hidden":false,"mode":"md"},{"id":543,"value":"

In this notebook, we evaluate the polynomial regression from data introduced by the user.

","pinned":false,"hidden":false,"mode":"html"},{"id":506,"value":"

For a given data set of x,y pairs, a polynomial regression of order n, f(x), is evaluated:

","pinned":false,"hidden":false,"mode":"html"},{"id":511,"value":"tex.block`\nf(x) = a_0 + a_1 \\, x + a_2 \\, x^2 + \\, a_3 \\, x^3 + \\cdots + a_n \\, x^n\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":524,"value":"where a₀, a₁, a₂, a₃, … , a_n, represent coefficients.","pinned":false,"hidden":false,"mode":"html"},{"id":366,"value":"

① Provide x,y data

\nEnter or paste free-form x,y paired data in the following window:^(*)

* Enter data: Use single space between x and y data, press enter to introduce a new point.\n
Paste data: Copy x,y data from a spreadsheet.\n

","pinned":false,"hidden":false,"mode":"html"},{"id":216,"value":"viewof data = html``","pinned":false,"hidden":false,"mode":"ojs"},{"id":561,"value":"
\n

② Define the polynome degree

","pinned":false,"hidden":false,"mode":"html"},{"id":189,"value":"viewof polynome_degree = Inputs.range([1, 10], {value: 4, step: 1, label: \"Polynome degree\"})","pinned":false,"hidden":false,"mode":"ojs"},{"id":576,"value":"
\n

③ Compare data with polynomial regression

","pinned":false,"hidden":false,"mode":"html"},{"id":115,"value":"Plot.plot({\n inset: 0,\n grid: true,\n height: 500,\n width: 700,\n x: {label: \"x → \"},\n y: {label: \"↑ y \"},\n style: { fontSize: \"11.5px\" },\n marks: [\n Plot.ruleX([Math.min(...x),Math.max(...x)]),\n Plot.ruleY([Math.min(...y),Math.max(...y)]),\n Plot.line(fitting_curve, {x: \"x\", y: \"y\", curve: \"catmull-rom\", stroke: color_stroke, strokeDasharray: \"7,4\"}),\n Plot.dot(data_curve, {x: \"x\", y: \"y\", curve: \"catmull-rom\", stroke: color_dot, fill: 'white', r: 4}), \n ]\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":990,"value":"viewof color_dot = Inputs.color({ label: \"Dot color\", value: \"#4472C4\", datalist: [\"#4472C4\", \"#176A00\", \"#FF8000\", \"#C50707\", \"#000000\"] })","pinned":false,"hidden":false,"mode":"ojs"},{"id":144,"value":"viewof color_stroke = Inputs.color({ label: \"Stroke color\", value: \"#FF8000\", datalist: [\"#FF8000\", \"#C50707\", \"#00B4FF\", \"#9C9C9C\", \"#000000\"] })","pinned":false,"hidden":false,"mode":"ojs"},{"id":766,"value":"
\n

Polynomial regression function:

","pinned":false,"hidden":false,"mode":"html"},{"id":480,"value":"$31","pinned":false,"hidden":false,"mode":"ojs"},{"id":814,"value":"

For calculation purposes, copy the following line:

","pinned":false,"hidden":false,"mode":"html"},{"id":777,"value":"{\n var polynome = \"y = \" + polynome_terms[0];\n \n for(var i=1; i<=polynome_degree; i++)\n {\n polynome_terms[i] >=0 ? polynome += \" + \" : polynome += \" - \";\n polynome += Math.abs(polynome_terms[i]) + \"*x\"\n i>=2 ? polynome += \"^\" +i : \"\";\n }\n return html`${polynome}`;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":65,"value":"---\n## Appendix","pinned":false,"hidden":false,"mode":"md"},{"id":212,"value":"### Data and parameters","pinned":false,"hidden":false,"mode":"md"},{"id":274,"value":"// x = [ 0, 10, 20,30, 40, 50, 60,70, 80, 90, 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200];\nx = {\n var data_x = [];\n var data_ = data.replaceAll(\" \",\"\") // data manipulation to make this notebook compatible\n // with data copied from https://observablehq.com/@christophe-yamahata/on-line-tool-calque\n data_ = data_.replaceAll(\" \",\"\\t\").split(\"\\n\");\n const n = data_.length;\n for(var i=0; i < n; i++){\n if(!(isNaN(parseFloat(data_[i].split('\\t')[0])) || isNaN(parseFloat(data_[i].split('\\t')[1])))) {\n data_x.push(parseFloat(data_[i].split('\\t')[0])); \n }\n }\n return data_x;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":52,"value":"// y = [ 0, 0.0003, 0.0015, 0.0053, 0.017, 0.045, 0.1, 0.21, 0.37, 0.6, 0.9, 1.2, 1.7, 2.1, 2.6, 3.2, 3.7, 4.3, 4.9, 5.5, 6.1 ]\ny = {\n var data_y = [];\n var data_ = data.replaceAll(\" \",\"\") // data manipulation to make this notebook compatible\n // with data copied from https://observablehq.com/@christophe-yamahata/on-line-tool-calque\n data_ = data_.replaceAll(\" \",\"\\t\").split(\"\\n\");\n const n = data_.length;\n for(var i=0; i < n; i++){\n if(!(isNaN(parseFloat(data_[i].split('\\t')[0])) || isNaN(parseFloat(data_[i].split('\\t')[1])))) {\n data_y.push(parseFloat(data_[i].split('\\t')[1])); \n }\n }\n return data_y;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":153,"value":"data_curve = {\n var data_curve = [];\n const n = x.length;\n for(var i = 0; i < n; i++) \n {\n data_curve.push({\"x\":x[i],\"y\":y[i]}); \n }\n return data_curve;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":137,"value":"fitting_curve = {\n var fitting_curve = [];\n var x_ = [];\n var y_ = [];\n const n = 500;\n const dx = (Math.max(...x)-Math.min(...x))/n;\n for(var i = 0; i < n; i++) \n {\n x_[i] = Math.min(...x) + dx*i;\n y_[i] = solver(x_[i]);\n fitting_curve.push({\"x\":x_[i],\"y\":y_[i]});\n }\n return fitting_curve;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":475,"value":"polynome_terms = poly.computeCoefficients(polynome_degree);","pinned":false,"hidden":false,"mode":"ojs"},{"id":72,"value":"### Javascript functions","pinned":false,"hidden":false,"mode":"md"},{"id":14,"value":"poly = Polyfit(x, y)","pinned":true,"hidden":false,"mode":"ojs"},{"id":33,"value":"solver = poly.getPolynomial(polynome_degree);","pinned":true,"hidden":false,"mode":"ojs"},{"id":26,"value":"$32","pinned":true,"hidden":false,"mode":"ojs"}],"files":"$W33","stdlib":"1","titling":"auto"}