31:T944,{ let dotplot = Plot.plot({ height: 500, width: 600, grid: true, marginLeft: 50, marginRight: 30, marginTop: 40, marginBottom: 70, insetTop: 0, insetBottom: 0, style: { fontSize: "13px", color: "#000" }, x: { label: " x → ", tickPadding: 10, labelOffset: 60, tickFormat: "f", tickSize: 9, ticks: 10, domain: [xMin,xMax], }, y: { label: "↑ y", tickPadding: 6, labelOffset: 30, tickFormat: "f", ticks: 10, tickSize: 6, domain: [yMin, yMax] }, caption: html`
Figure: Comparison of data points with joint optimization results.
⚪ Data points
– – – Joint optimization without bounds
—— Joint optimization with bounds
`, marks: [ Plot.frame({stroke: "black"}), Plot.ruleY([0]), Plot.dot(dataPoints.filter(point => point.dataset === 'y1'), {x: "x", y: "y", stroke:"#111", fill: strokeColor1, strokeOpacity: 0.4, r: dotSize}), Plot.dot(dataPoints.filter(point => point.dataset === 'y2'), {x: "x", y: "y", stroke:"#111", fill: strokeColor2, strokeOpacity: 0.4, r: dotSize}), Plot.line(createCurve(optimizedParamsNoConstraints,'y1'), {x: "x", y: "y", curve: "catmull-rom", stroke: strokeColor1, strokeWidth: 0.8, strokeOpacity: 0.8, strokeDasharray: "4,4"}), Plot.line(createCurve(optimizedParamsNoConstraints,'y2'), {x: "x", y: "y", curve: "catmull-rom", stroke: strokeColor2, strokeWidth: 0.8, strokeOpacity: 0.8, strokeDasharray: "4,4"}), Plot.line(createCurve(optimizedParamsBounded,'y1'), {x: "x", y: "y", curve: "catmull-rom", stroke: strokeColor1, strokeWidth: 0.8, strokeOpacity: 0.8}), Plot.line(createCurve(optimizedParamsBounded,'y2'), {x: "x", y: "y", curve: "catmull-rom", stroke: strokeColor2, strokeWidth: 0.8, strokeOpacity: 0.8}), ] }) // Label in bold/large font let select_labels = d3.select(dotplot) .selectAll("svg > g > text") .style("font-weight", "bold") .style("font-size", "16px"); return dotplot; }32:T964,// Modified fminsearch function with bounds: Please note that this implementation is not optimized. // Adapted from: https://github.com/jonasalmeida/fminsearch // See associated notebook: https://observablehq.com/@christophe-yamahata/nonlinear-regression-with-bounds-in-javascript 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 if (!Opt.bounds) { Opt.bounds = Parm0.map(function() { return [-Infinity, Infinity] }) }; // Default no bounds 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; let step = cloneVector(Opt.step); // Ensure step is reinitialized var funParm = function(P) { return Opt.objFun(y, fun(x, P)) } //function (of Parameters) to minimize // Apply bounds var applyBounds = function(P, bounds) { return P.map(function(p, i) { return Math.max(bounds[i][0], Math.min(bounds[i][1], p)); }); } // Silly multi-univariate screening for (var i = 0; i < Opt.maxIter; i++) { for (var j = 0; j < n; j++) { // take a step for each parameter P1 = cloneVector(P0); P1[j] += step[j]; P1 = applyBounds(P1, Opt.bounds); // Apply bounds after updating the parameter if (funParm(P1) < funParm(P0)) { // if parm value going in the right direction step[j] = 1.2 * step[j]; // then go a little faster P0 = cloneVector(P1); } else { step[j] = -(0.5 * step[j]); // otherwise reverse and go slower } } if (Opt.display) { if (i > (Opt.maxIter - 10)) { console.log(i + 1, funParm(P0), P0) } } } return P0 }33:T7df,python_script = ` ''' This Python script demonstrates joint optimization using least squares regression with bounded parameters. The goal is to find the optimal parameters for related datasets by minimizing the sum of squared differences between observed and predicted values, while ensuring that the parameters stay within specified bounds. Source with details about the example: https://observablehq.com/@christophe-yamahata/joint-optimization-with-bounded-parameters Author: C. Yamahata Date: May 2024 ''' import numpy as np from scipy.optimize import minimize # Define the model function def model(parameters, x, dataset): a, b, c = parameters # Compute predicted values based on dataset labels return np.where(dataset == 'y1', a * x + b, a * x + c) # Define the objective function for least squares def objective(parameters, x, y, dataset): y_pred = model(parameters, x, dataset) # Return the sum of squared differences return np.sum((y - y_pred) ** 2) # Define the data for dataset 1 x_1 = np.array([0.1, 2.1, 3.15, 4.05, 5.1, 6.14, 7.3, 8.25, 9.15, 10.1]) y_1 = np.array([5.485, 11.486, 15.764, 22.448, 26.204, 32.339, 37.523, 43.72, 47.206, 51.495]) # Define the data for dataset 2 x_2 = np.array([0.12, 1.12, 2.22, 3.4, 4.1, 5.05, 6.02, 7.12, 8.1, 9.17, 10.05]) y_2 = np.array([13.674, 19.861, 30.347, 32.001, 40.433, 38.473, 46.0411, 55.523, 55.066, 65.012, 66.704]) # Concatenate the datasets x_all = np.concatenate((x_1, x_2)) y_all = np.concatenate((y_1, y_2)) dataset_all = np.array(['y1'] * len(x_1) + ['y2'] * len(x_2)) # Initial guess for parameters [a, b, c] initial_parameters = [5.2, 1.3, 12] # Define bounds for parameters a, b, and c bounds = [(5.1948, 5.2052), (1.04, 1.56), (10.8, 13.2)] # Perform optimization with bounds optimized_params_bounded = minimize(objective, initial_parameters, args=(x_all, y_all, dataset_all), bounds=bounds) # Output the optimized parameters print("Optimized Parameters with Bounds:\\n", optimized_params_bounded.x) `34:[{"id":"49c05b3e889c7e17e2c0d45a5788f60eb69fe3d8790210c693825d33a2030c7ed6da0b92109b0047df1c16cdc5b3ce50ab6a9bbbd367a3fd1fe0552634a408dd","href":"https://static.observableusercontent.com/files/49c05b3e889c7e17e2c0d45a5788f60eb69fe3d8790210c693825d33a2030c7ed6da0b92109b0047df1c16cdc5b3ce50ab6a9bbbd367a3fd1fe0552634a408dd","name":"Python_output.PNG","type":"image/png","lastModified":1717176488169.751,"size":2033}] 2f:{"title":"Joint Optimization with Bounded Parameters","cells":[{"id":0,"value":"# Joint Optimization with Bounded Parameters","pinned":false,"hidden":false,"mode":"md"},{"id":922,"value":"

\n This notebook demonstrates joint optimization using least squares regression with bounded parameters. The goal is to find the optimal parameters for related datasets by minimizing the sum of squared differences between observed and predicted values, while ensuring that the parameters stay within specified bounds.\n

Joint Optimization – Simultaneously optimize parameters for several related datasets.
Least Squares Regression – Minimize the error between observed and predicted values using the least squares method.
Bounded Parameters – Impose limits on the parameters to ensure they remain within specified bounds during optimization.

\n\n

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

Example usage

\n For illustration purposes, consider two datasets defined by the following functions:\n
${tex`y_1(x) = a \\times x + b`}\n
${tex`y_2(x) = a \\times x + c`}\n

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

Here, ${tex`a`}, ${tex`b`}, and ${tex`c`}, are the fitting parameters to be determined.

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

Data points

\n The data points ${tex`(x_1,y_1)`} and ${tex`(x_2,y_2)`} are stored in separate arrays:\n

","pinned":false,"hidden":false,"mode":"html"},{"id":12,"value":"x_1 = [0.1,2.1,3.15,4.05,5.1,6.14,7.3,8.25,9.15,10.1]","pinned":false,"hidden":false,"mode":"ojs"},{"id":14,"value":"y_1 = [5.485, 11.486, 15.764, 22.448, 26.204, 32.339, 37.523, 43.720, 47.206, 51.495, 56.697 ]","pinned":false,"hidden":false,"mode":"ojs"},{"id":1648,"value":"x_2 = [0.12, 1.12, 2.22, 3.4, 4.1, 5.05, 6.02, 7.12, 8.1, 9.17, 10.05]","pinned":false,"hidden":false,"mode":"ojs"},{"id":1651,"value":"y_2 = [13.674, 19.861, 30.347, 32.001, 40.433, 38.473, 46.0411, 55.523, 55.066, 65.012, 66.704 ]","pinned":false,"hidden":false,"mode":"ojs"},{"id":1657,"value":"

To facilitate data manipulation, we consolidate the data points into a single array, dataPoints:

","pinned":false,"hidden":false,"mode":"html"},{"id":1654,"value":"// Combine the data points into a single array\ndataPoints = [\n ...x_1.map((x, i) => ({ x, y: y_1[i], dataset: 'y1' })),\n ...x_2.map((x, i) => ({ x, y: y_2[i], dataset: 'y2' }))\n]","pinned":true,"hidden":false,"mode":"ojs"},{"id":1846,"value":"
\n

Graphic illustration

\n In the following graph, you can compare the final results:\n

Joint optimization without bounding limits ( - - - )
Joint optimization with bounding limits ( — )

\n
","pinned":false,"hidden":false,"mode":"html"},{"id":1739,"value":"$31","pinned":false,"hidden":false,"mode":"ojs"},{"id":1853,"value":"

","pinned":false,"hidden":false,"mode":"html"},{"id":1664,"value":"
\n

Model

From the equations ${tex`y_1(x)`} and ${tex`y_2(x)`} given above, we define the model as follows:

","pinned":false,"hidden":false,"mode":"html"},{"id":18,"value":"// Define the model function\nfunction model(params, x, dataset) {\n const [a, b, c] = params;\n if (dataset === 'y1') {\n return a * x + b;\n } else if (dataset === 'y2') {\n return a * x + c;\n }\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":1689,"value":"

Initial guess bounding limits, number of iterations

The initial guess parameters are loaded into the array `initialParams`:

","pinned":false,"hidden":false,"mode":"html"},{"id":1687,"value":"// Initial guess for parameters [a, b, c]\ninitialParams = [5.2, 1.3, 12]","pinned":false,"hidden":false,"mode":"ojs"},{"id":1702,"value":"

Lower and upper bounds for parameters are stored in the array `bounds`.

","pinned":false,"hidden":false,"mode":"html"},{"id":1693,"value":"// Define bounds for parameters a, b, and c\nbounds = {\n return [\n [0.999 * initialParams[0], 1.001 * initialParams[0]], // Bounds for a \n [0.8 * initialParams[1], 1.2 * initialParams[1]], // Bounds for b\n [0.9 * initialParams[2], 1.1 * initialParams[2]], // Bounds for c\n ];\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":1750,"value":"

Value for maximum number of iterations is stored in the variable `maxIterations`:

","pinned":false,"hidden":false,"mode":"html"},{"id":1712,"value":"maxIterations = 200","pinned":false,"hidden":false,"mode":"ojs"},{"id":1780,"value":"

Parameters `bounds` and `maxIterations` are stored in the object `Opt` which is used by the optimization function fminsearch.

","pinned":false,"hidden":false,"mode":"html"},{"id":1715,"value":"Opt = {\n return {\n maxIter: maxIterations, // Optional: Specify the maximum number of iterations\n bounds: bounds, // Bounds for the parameters\n }\n};","pinned":false,"hidden":false,"mode":"ojs"},{"id":1777,"value":"
\n

Calculation of the parameters with bounding limits

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

\n The least squares regression with bounding limits is completed using the JavaScript function fminsearch as follows:\n

\noptimizedParamsBounded = fminsearch(model,initParams,x,y, Opt)\n","pinned":false,"hidden":false,"mode":"html"},{"id":1791,"value":"

The resulting optimized parameters with bounding limits are:

\n ${tex`a = `} ${optimizedParamsBounded[0].toFixed(4)}\n
${tex`b = `} ${optimizedParamsBounded[1].toFixed(4)}\n
${tex`c = `} ${optimizedParamsBounded[2].toFixed(4)}","pinned":false,"hidden":false,"mode":"html"},{"id":1795,"value":"optimizedParamsBounded = fminsearch((x, params) => dataPoints.map(p => model(params, p.x, p.dataset)), initialParams, dataPoints.map(p => p.x), dataPoints.map(p => p.y),Opt);","pinned":false,"hidden":false,"mode":"ojs"},{"id":845,"value":"
\n

Calculation of the parameters without bounding limits

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

\n This result can be compared to the solution obtained when the function fminsearch is called without the optional bounds:\n

\noptimizedParamsNoConstraints = fminsearch(model,initParams,x,y, {maxIter: maxIterations})\n","pinned":false,"hidden":false,"mode":"html"},{"id":1525,"value":"

The resulting optimized parameters without bounding limits are:

\n ${tex`a = `} ${optimizedParamsNoConstraints[0].toFixed(4)}\n
${tex`b = `} ${optimizedParamsNoConstraints[1].toFixed(4)}\n
${tex`c = `} ${optimizedParamsNoConstraints[2].toFixed(4)}","pinned":false,"hidden":false,"mode":"html"},{"id":1742,"value":"optimizedParamsNoConstraints = fminsearch((x, params) => dataPoints.map(p => model(params, p.x, p.dataset)), initialParams, dataPoints.map(p => p.x), dataPoints.map(p => p.y));\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":2039,"value":"

Equivalent Python script

\nThe following script provides a similar implementation in Python. The script uses the SciPy library.\n

","pinned":false,"hidden":false,"mode":"html"},{"id":2042,"value":"Inputs.button(\"Copy Python script to clipboard\", {value: null, reduce: () => navigator.clipboard.writeText(python_script)})","pinned":false,"hidden":false,"mode":"ojs"},{"id":2061,"value":"
\n

\n
Results generated by the Python script.\n