31:T654,

Very frequently, one wishes to fit some experimental data to a simple nonlinear analytical model. This can be completed through nonlinear regression. This notebook shows how this can be implemented using Javascript. Furthermore, in the files pane , several examples can be downloaded and tested in GNU Octave.

32:T7be,

Nonlinear regression in Octave

A nonlinear fitting algorithm can be computed in Matlab or in GNU Octave using one of the following functions: fminunc, fminsearch.
You can click on the files pane icon on the right-hand side or click on the following links to download the

CSV data file: "UV-LED_spectrum.csv"
Octave script: "nonlinear_fitting_example.m"

Note that you will need to place the *.csv file in the same folder as the *.m file to execute the script in GNU Octave.

33:T5e3,You can click on the files pane icon on the right-hand side or click on the following links to download the

CSV data file: "UV-LED_spectrum.csv"
Python script: "nonlinear_fitting_example.py"

34:T56b,// Source: https://github.com/jonasalmeida/fminsearch fminsearch = function(fun,Parm0,x,y,Opt){ if(!Opt){Opt={}}; if(!Opt.maxIter){Opt.maxIter=1000}; if(!Opt.step){// initial step is 1/100 of initial value (remember not to use zero in Parm0) Opt.step=Parm0.map(function(p){return p/100}); Opt.step=Opt.step.map(function(si){if(si==0){return 1}else{ return si}}); // convert null steps into 1's }; if(typeof(Opt.display)=='undefined'){Opt.display=true}; if(!Opt.objFun){Opt.objFun=function(y,yp){return y.map(function(yi,i){return Math.pow((yi-yp[i]),2)}).reduce(function(a,b){return a+b})}} //SSD var cloneVector=function(V){return V.map(function(v){return v})}; var ya,y0,yb,fP0,fP1; var P0=cloneVector(Parm0),P1=cloneVector(Parm0); var n = P0.length; var step=Opt.step; var funParm=function(P){return Opt.objFun(y,fun(x,P))}//function (of Parameters) to minimize // silly multi-univariate screening for(var i=0;i(Opt.maxIter-10)){console.log(i+1,funParm(P0),P0)}} } return P0 }35:T89f,python_script = ` """ Nonlinear Regression Example in Python This script performs nonlinear regression using the SciPy library to fit a model to observed data points. For further details and a similar implementation in JavaScript, visit: https://observablehq.com/@christophe-yamahata/nonlinear-regression-in-javascript You can download the accompanying CSV file "UV-LED_spectrum.csv" from the above URL to use as input data. Author: C. Yamahata Date: May 2024 """ import numpy as np import pandas as pd from scipy.optimize import least_squares # Define the model function def model(x, p): return p[0] * np.exp(-0.5 * ((x - p[1]) / p[2])**2) + p[3] # Uncomment the following lines to test with CSV data # return p[0]* np.exp(-0.5*np.power((x-p[1])/p[2],2)) \\ # + p[3]*np.exp(-0.5*np.power((x-p[4])/p[5],2)) \\ # + p[6]*np.exp(-0.5*np.power((x-p[7])/p[8],2)) + p[9]*np.exp(-0.5*np.power((x-p[10])/p[11],2)) # Define the objective function def objective_function(p, x, y): return model(x, p) - y # Given data x = np.array([32,37,42,47,52,57,62,67,72,77,82,87,92]) y = np.array([749,1525,1947,2201,2380,2537,2671,2758,2803,2943,3007,2979,2992]) p_init = np.array([-70000,-90,40,3000]) # Uncomment the following lines to test with CSV data # df = pd.read_csv("UV-LED_spectrum.csv", delimiter=',', skiprows=1) # Read data from CSV file # x = df.iloc[:, 0].values # Extract x and y values from DataFrame # y = df.iloc[:, 1].values # p_init = np.array([0.9, 365, 5, 0.2, 390, 15, 0.4, 400, 10, 0.2, 410, 15]) # Initial guess # Define lower and upper bounds for each parameter lower_bounds = [-np.inf, -np.inf, -np.inf, -np.inf] upper_bounds = [np.inf, np.inf, np.inf, np.inf] # Uncomment the following lines to test with CSV data # lower_bounds = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0] # upper_bounds = [np.inf, np.inf, np.inf, np.inf, np.inf, np.inf, np.inf, np.inf, np.inf, np.inf, np.inf, np.inf] bounds = (lower_bounds, upper_bounds) # Optimize parameters using least squares with bounds result = least_squares(objective_function, p_init, bounds=bounds, args=(x, y)) print("\\nBest-fit parameters:") for parameter in result.x: print(f"{parameter}") `36:[{"id":"43b94783684dc404cdfd5003206d23df73207827af027586c1d95b0755965479f4a961f433d3ff315cf801a228af3c0b10a86ab01a9506cc1bcd3cab9b570006","href":"https://static.observableusercontent.com/files/43b94783684dc404cdfd5003206d23df73207827af027586c1d95b0755965479f4a961f433d3ff315cf801a228af3c0b10a86ab01a9506cc1bcd3cab9b570006","name":"nonlinear_fitting_example.py","type":"application/octet-stream","lastModified":1716644024553.62,"size":2239},{"id":"39cfe0b41c3af579aa5067a9602a5f6df43291e742c21e9d4f6a86e6a4728ff40d1acbdff561fdc800fd93cd2d0c08fb2e635b402bb344fe8365fd5948c5f96b","href":"https://static.observableusercontent.com/files/39cfe0b41c3af579aa5067a9602a5f6df43291e742c21e9d4f6a86e6a4728ff40d1acbdff561fdc800fd93cd2d0c08fb2e635b402bb344fe8365fd5948c5f96b","name":"nonlinear_fitting_example.m","type":"application/octet-stream","lastModified":1667750952888.179,"size":1503},{"id":"77963be9a6117a0fff76b15d9f5c10d4c61c3a46c9f43c49ea9b12ce25a2c0912c5dcfdf1ce6bcecec2c33fa48fcb9ee6d5365ba02461c6211c6b0efee575b89","href":"https://static.observableusercontent.com/files/77963be9a6117a0fff76b15d9f5c10d4c61c3a46c9f43c49ea9b12ce25a2c0912c5dcfdf1ce6bcecec2c33fa48fcb9ee6d5365ba02461c6211c6b0efee575b89","name":"UV-LED_spectrum.csv","type":"text/csv","lastModified":1667750112981.163,"size":4942},{"id":"2a399078b0cf6910eeec87b38887c24340bd37770a06fa08fad3518aa3e6059e460d10955751d6c714a9684f82eee045f5c2bba74cc06cc1ab971394f5f88608","href":"https://static.observableusercontent.com/files/2a399078b0cf6910eeec87b38887c24340bd37770a06fa08fad3518aa3e6059e460d10955751d6c714a9684f82eee045f5c2bba74cc06cc1ab971394f5f88608","name":"nonlinear_fitting_example_8.m","type":"application/octet-stream","lastModified":1668029946973.408,"size":3017},{"id":"62df340aedf06e273c8707ca73b535b7e335202cd5a96410f10ddfcefd0892d6122c25be2f204209aabf06158a3728717a12571b0414da7f89f0db1ac0a2950e","href":"https://static.observableusercontent.com/files/62df340aedf06e273c8707ca73b535b7e335202cd5a96410f10ddfcefd0892d6122c25be2f204209aabf06158a3728717a12571b0414da7f89f0db1ac0a2950e","name":"nonlinear_fitting_example_9.m","type":"application/octet-stream","lastModified":1668032794253.18,"size":1183}] 2f:{"title":"Nonlinear regression in JavaScript","cells":[{"id":0,"value":"# Nonlinear regression in JavaScript","pinned":false,"hidden":false,"mode":"md"},{"id":922,"value":"$31","pinned":false,"hidden":false,"mode":"html"},{"id":558,"value":"introduction = md`

\n As an illustrative example, we will consider the spectral power distribution of a typical UV-LED with two emission peaks (365 nm and 405 nm in the example given hereafter).\n The spectral power distribution can be approximated as a sum of Gaussian functions: ^{[[1](#Yamahata2020)]}

\n `","pinned":false,"hidden":false,"mode":"ojs"},{"id":37,"value":"tex.block`\n U(\\lambda) = \\sum_{i=0}^{\\infty} A_i e^{-\\frac{1}{2} (\\frac{\\lambda-\\lambda_i}{\\sigma_i})^2 }\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":642,"value":"

\n Fortunately, an excellent approximation can already be obtained by summing only a limited number of Gaussian curves. Hereafter, we will use only three Gaussian curves: ${tex`i=0,1,2,3`}.\n

","pinned":false,"hidden":false,"mode":"html"},{"id":693,"value":"$32","pinned":false,"hidden":false,"mode":"html"},{"id":1116,"value":"
\n

Nonlinear regression in Python

\nSimilarly, the following Python script performs nonlinear regression using the SciPy library to fit a model to observed data points.","pinned":false,"hidden":false,"mode":"html"},{"id":1106,"value":"Inputs.button(\"Copy Python script to clipboard\", {value: null, reduce: () => navigator.clipboard.writeText(python_script)})","pinned":false,"hidden":false,"mode":"ojs"},{"id":1087,"value":"$33","pinned":false,"hidden":false,"mode":"html"},{"id":665,"value":"
\n

Nonlinear regression in JavasScript

In this notebook, we will test a similar algorithm that can be directly executed in a web browser:\n @jonasalmeida has implemented a Javascript function (fminsearch.js) that can be found on GitHub.\n

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

An example of use is shown below:

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

\n
\n

\n    x = [32,37,42,47,52,57,62,67,72,77,82,87,92]    
\n    y = [749,1525,1947,2201,2380,2537,2671,2758,2803,2943,3007,2979,2992]   
\n    model = function(x,p){return x.map(function(x_i){return ( 
\n       p[0]*Math.exp(-0.5*(Math.pow((x_i-p[1])/p[2],2))) + p[3]    
 \n    )})}  
\n    p_init = [-70000,-90,40,3000];  
\n    fitting = fminsearch(model,p_init,x,y)\n

\n

\n

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

Let's now adapt this code to find the fitting parameters in the case of the spectral power distribution of a UV-LED whose measurement data are given in this CSV file.

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

\n First, the x (wavelength, ${tex`\\lambda`}) and y (amplitude, ${tex`U`}) data are retrieved from the CSV file and stored in two arrays:\n

","pinned":false,"hidden":false,"mode":"html"},{"id":12,"value":"// x = [32,37,42,47,52,57,62,67,72,77,82,87,92]\n\n// data: x = lambda\nx = {\n let lambda = [];\n for (let i = 0; i < LED_spectrum.length; i++) {\n lambda.push(LED_spectrum[i].lambda);\n }\n return lambda;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":14,"value":"// y = [749,1525,1947,2201,2380,2537,2671,2758,2803,2943,3007,2979,2992]\n\n// data: y = U\n y = {\n let U = [];\n for (let i = 0; i < LED_spectrum.length; i++) {\n U.push(LED_spectrum[i].U);\n }\n return U;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":886,"value":"\n

Then, one needs to define the model. From the analytical model given above, we write:

","pinned":false,"hidden":false,"mode":"html"},{"id":18,"value":"// model = function(x,p){return x.map(function(x_i){return (\n// p[0]*Math.exp(-0.5*(Math.pow((x_i-p[1])/p[2],2))) + p[3] \n// )})}\nmodel = function(x,p){return x.map(function(x_i){return (\n p[0]*Math.exp(-0.5*Math.pow((x_i-p[1])/p[2],2)) + \n p[3]*Math.exp(-0.5*Math.pow((x_i-p[4])/p[5],2)) +\n p[6]*Math.exp(-0.5*Math.pow((x_i-p[7])/p[8],2)) + \n p[9]*Math.exp(-0.5*Math.pow((x_i-p[10])/p[11],2))\n)})}","pinned":true,"hidden":false,"mode":"ojs"},{"id":786,"value":"

Next, the initialisation parameters (i.e., initial guess) are loaded in the array p_init:

","pinned":false,"hidden":false,"mode":"html"},{"id":28,"value":"// p_init = [-70000,-90,40,3000];\np_init = [0.9, 365, 5, 0.2, 390, 15, 0.4, 400, 10, 0.2, 410, 15]; ","pinned":true,"hidden":false,"mode":"ojs"},{"id":790,"value":"

\n Finally, the nonlinear regression is completed using fminsearch:\n

","pinned":false,"hidden":false,"mode":"html"},{"id":22,"value":"fitting = fminsearch(model,p_init,x,y)","pinned":true,"hidden":false,"mode":"ojs"},{"id":901,"value":"
\n

In the following graph, the final result can be compared to the initial estimate.

","pinned":false,"hidden":false,"mode":"html"},{"id":61,"value":"Plot.plot({\n inset: 0,\n grid: true,\n height: 450,\n width: 800,\n x: {label: \"Wavelength, λ (nm) → \"},\n y: {label: \"↑ y (normalized data)\"},\n style: { fontSize: \"11.5px\" },\n marks: [\n Plot.ruleX([x_min,x_max]),\n Plot.ruleY([y_min,y_max]),\n \n Plot.line(initial_curve, {x: \"x\", y: \"y\", curve: \"catmull-rom\", stroke: color_initial, strokeOpacity: 0.3}), \n \n Plot.line(fitting_curve, {x: \"x\", y: \"y\", curve: \"catmull-rom\", stroke: color_stroke}),\n \n Plot.dot(LED_spectrum, {x: \"lambda\", y: \"U\", curve: \"catmull-rom\", fill: \"white\", stroke: color_dot, r: 1.6, strokeOpacity: 0.7}), \n\n\n ]\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":564,"value":"---\n## Reference","pinned":false,"hidden":false,"mode":"md"},{"id":567,"value":"Yamahata2020 = reference(\"introduction\", 1, \"C. Yamahata, \\\"Broadband radiometric measurement of light-emitting diodes with a photodiode,\\\" White Paper released by idonus Sàrl, December 2020.
DOI: 10.13140/RG.2.2.34765.05600/2\")","pinned":false,"hidden":false,"mode":"ojs"},{"id":510,"value":"---\n## Nonlinear regression: JavaScript function","pinned":false,"hidden":false,"mode":"md"},{"id":5,"value":"$34","pinned":true,"hidden":false,"mode":"ojs"},{"id":25,"value":"---\n## Appendix","pinned":false,"hidden":false,"mode":"md"},{"id":656,"value":"url_octave_example = await FileAttachment(\"nonlinear_fitting_example.m\").url();","pinned":false,"hidden":false,"mode":"ojs"},{"id":1091,"value":"url_python_example = await FileAttachment(\"nonlinear_fitting_example.py\").url();","pinned":false,"hidden":false,"mode":"ojs"},{"id":1109,"value":"$35","pinned":false,"hidden":false,"mode":"ojs"},{"id":700,"value":"url_csv_data = await FileAttachment(\"UV-LED_spectrum.csv\").url();","pinned":false,"hidden":false,"mode":"ojs"},{"id":308,"value":"LED_spectrum = FileAttachment(\"UV-LED_spectrum.csv\").csv({typed: true})","pinned":false,"hidden":false,"mode":"ojs"},{"id":157,"value":"initial_curve = {\n var initial_curve = [];\n var x_ = [];\n var y_ = [];\n const n = 500;\n const dx = (x_max-x_min)/n;\n for(var i = 0; i < n; i++) \n {\n x_[i] = x_min + dx*i;\n y_[i] = p_init[0]*Math.exp(-0.5*Math.pow((x_[i]-p_init[1]) / p_init[2],2)) + p_init[3]*Math.exp(-0.5*Math.pow((x_[i]-p_init[4]) / p_init[5],2)) +\n p_init[6]*Math.exp(-0.5*Math.pow((x_[i]-p_init[7]) / p_init[8],2)) + p_init[9]*Math.exp(-0.5*Math.pow((x_[i]-p_init[10]) / p_init[11],2))\n initial_curve.push({\"x\":x_[i],\"y\":y_[i]});\n }\n\n return initial_curve;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":134,"value":"fitting_curve = {\n var fitting_curve = [];\n var x_ = [];\n var y_ = [];\n const n = 500;\n const dx = (x_max-x_min)/n;\n for(var i = 0; i < n; i++) \n {\n x_[i] = x_min + dx*i;\n y_[i] = fitting[0]*Math.exp(-0.5*Math.pow((x_[i]-fitting[1])/fitting[2],2))\n + fitting[3]*Math.exp(-0.5*Math.pow((x_[i]-fitting[4])/fitting[5],2)) \n + fitting[6]*Math.exp(-0.5*Math.pow((x_[i]-fitting[7])/fitting[8],2)) \n + fitting[9]*Math.exp(-0.5*Math.pow((x_[i]-fitting[10])/fitting[11],2))\n\n fitting_curve.push({\"x\":x_[i],\"y\":y_[i]});\n }\n\n return fitting_curve;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":122,"value":"color_dot = \"#4472C4\"","pinned":false,"hidden":false,"mode":"ojs"},{"id":454,"value":"color_stroke = \"orange\"","pinned":false,"hidden":false,"mode":"ojs"},{"id":465,"value":"color_initial = \"grey\"","pinned":false,"hidden":false,"mode":"ojs"},{"id":287,"value":"x_min = 340;","pinned":false,"hidden":false,"mode":"ojs"},{"id":291,"value":"x_max = 450;","pinned":false,"hidden":false,"mode":"ojs"},{"id":339,"value":"y_min = 0;","pinned":false,"hidden":false,"mode":"ojs"},{"id":341,"value":"y_max = 1.02;","pinned":false,"hidden":false,"mode":"ojs"},{"id":569,"value":"// https://observablehq.com/@jashkenas/markdown-style-footnotes\nreference = (referencingCellName, number, text) => {\n return html`

↫ ${text}

`;\n}","pinned":false,"hidden":false,"mode":"ojs"}],"files":"$W36","stdlib":"1","titling":"auto"}