MISC_PIECEWISE

Overview

The MISC_PIECEWISE function fits a collection of specialized piecewise and physics-inspired models to data using nonlinear least squares optimization. This function leverages SciPy’s curve_fit to estimate model parameters that minimize the sum of squared residuals between the fitted function and observed data.

Piecewise linear models are particularly useful when data exhibits distinct behavioral regimes separated by breakpoints. The function supports two-segment and three-segment linear models with automatic breakpoint detection, where the transition point (or points) between segments is fitted along with the slopes and intercepts. The Heaviside Step Two Level model captures abrupt transitions between two constant values, useful for threshold detection in experimental data.

The function also includes physics-based models from magnetism theory. The Langevin Paramagnetic and Langevin Scaled Field models implement the Langevin function, which describes how magnetic dipoles align under an applied field in classical paramagnetic materials. Named after French physicist Paul Langevin who developed the theory in 1905, the function is defined as:

L(x) = \coth(x) - \frac{1}{x}

This function smoothly transitions from linear behavior at small fields to saturation at large fields, making it valuable for modeling magnetic hysteresis and saturation effects.

For multivariate outputs, the function provides the Planar Surface Two Outputs model (fitting two linear functions simultaneously) and the Linear Exponential Multivariate model (combining linear and exponential components across two output dimensions). These are useful when the dependent variable has multiple components that share the same independent variable.

The optimization uses the Levenberg-Marquardt algorithm by default, with automatic initial parameter estimation to improve convergence. The function returns fitted parameter values along with standard errors derived from the covariance matrix when available.

This example function is provided as-is without any representation of accuracy.

Excel Usage

=MISC_PIECEWISE(xdata, ydata, piecewise_model)

• xdata (list[list], required): The xdata value
• ydata (list[list], required): The ydata value
• piecewise_model (str, required): The piecewise_model value

Returns (list[list]): 2D list [param_names, fitted_values, std_errors], or error string.

Examples

Example 1: Demo case 1

Inputs:

piecewise_model	xdata	ydata
two_segment_linear_breakpoint	0	1.5
	1.5	4.8
	3	8.1
	3.5	9.2
	4.5	10
	5.5	10.8
	6.5	11.6

Excel formula:

=MISC_PIECEWISE("two_segment_linear_breakpoint", {0;1.5;3;3.5;4.5;5.5;6.5}, {1.5;4.8;8.1;9.2;10;10.8;11.6})

Expected output:

a1	k1	xi	k2
1.5	2.2	3.5	0.8
2.083e-15	8.593e-16	2.854e-15	1.664e-15

Example 2: Demo case 2

Inputs:

piecewise_model	xdata	ydata
three_segment_linear_breakpoint	0	0.5
	1.2	1.7
	2.4	2.9
	3.6	3.44
	4.8	3.92
	6	4.95
	7.2	6.75
	8.4	8.55

Excel formula:

=MISC_PIECEWISE("three_segment_linear_breakpoint", {0;1.2;2.4;3.6;4.8;6;7.2;8.4}, {0.5;1.7;2.9;3.44;3.92;4.95;6.75;8.55})

Expected output:

a1	k1	xi1	k2	xi2	k3
0.5	1	2.5	0.4	5.5	1.5
0	0	0	0	0	0

Example 3: Demo case 3

Inputs:

piecewise_model	xdata	ydata
heaviside_step_two_level	0	2.5
	1	2.5
	2	2.5
	3	5
	4	5
	5	5

Excel formula:

=MISC_PIECEWISE("heaviside_step_two_level", {0;1;2;3;4;5}, {2.5;2.5;2.5;5;5;5})

Expected output:

A	B	x1
2.5	5	2.5

Example 4: Demo case 4

Inputs:

piecewise_model	xdata	ydata
langevin_paramagnetic	0.2	-0.760716189912
	1.1	-0.19025145403
	2	0.460698113858
	3.5	1.282091986942
	5	1.675444725654
	6.5	1.86019989127

Excel formula:

=MISC_PIECEWISE("langevin_paramagnetic", {0.2;1.1;2;3.5;5;6.5}, {-0.760716189912;-0.19025145403;0.460698113858;1.282091986942;1.675444725654;1.86019989127})

Expected output:

y0	xc	C
0.1	1.5	2.2
7.358e-7	0.000001264	6.928e-7

Example 5: Demo case 5

Inputs:

piecewise_model	xdata	ydata
langevin_scaled_field	0.2	-0.140701683384
	1.4	0.145836567552
	2.6	0.423906840147
	3.8	0.67176890878
	5	0.877979751168
	6.2	1.041693485784

Excel formula:

=MISC_PIECEWISE("langevin_scaled_field", {0.2;1.4;2.6;3.8;5;6.2}, {-0.140701683384;0.145836567552;0.423906840147;0.67176890878;0.877979751168;1.041693485784})

Expected output:

y0	xc	C	s
0.05	1	1.8	2.5
0.00000259	0.00001105	0.000005469	0.000006939

Example 6: Demo case 6

Inputs:

piecewise_model	xdata	ydata
planar_surface_two_outputs	0	1	2
	1	1.75	1.5
	2	2.5	1
	3	3.25	0.5
	4	4	0
	5	4.75	-0.5

Excel formula:

=MISC_PIECEWISE("planar_surface_two_outputs", {0;1;2;3;4;5}, {1,2;1.75,1.5;2.5,1;3.25,0.5;4,0;4.75,-0.5})

Expected output:

a	b	c	d
1	0.75	2	-0.5
0	0	0	0

Example 7: Demo case 7

Inputs:

piecewise_model	xdata	ydata
linear_exponential_multivariate	0	0.5	2
	0.8	1.3	2.216999356
	1.6	2.1	2.492860338
	2.4	2.9	2.843547305
	3.2	3.7	3.289357351
	4	4.5	3.856094251

Excel formula:

=MISC_PIECEWISE("linear_exponential_multivariate", {0;0.8;1.6;2.4;3.2;4}, {0.5,2;1.3,2.216999356;2.1,2.492860338;2.9,2.843547305;3.7,3.289357351;4.5,3.856094251})

Expected output:

a1	a2	b1	b2	c
0.5	1.2	1	0.8	0.3
1.848e-7	0.000001529	7.628e-8	0.000001368	2.961e-7

Python Code

import numpy as np
from scipy.optimize import curve_fit as scipy_curve_fit
import math

def misc_piecewise(xdata, ydata, piecewise_model):
    """
    Fits misc_piecewise models to data using scipy.optimize.curve_fit. See https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html for details.

    See: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html

    This example function is provided as-is without any representation of accuracy.

    Args:
        xdata (list[list]): The xdata value
        ydata (list[list]): The ydata value
        piecewise_model (str): The piecewise_model value Valid options: Two Segment Linear Breakpoint, Three Segment Linear Breakpoint, Heaviside Step Two Level, Langevin Paramagnetic, Langevin Scaled Field, Planar Surface Two Outputs, Linear Exponential Multivariate.

    Returns:
        list[list]: 2D list [param_names, fitted_values, std_errors], or error string.
    """
    def _validate_data(xdata, ydata):
        """Validate and convert input data to numpy arrays."""

        def _convert_x(arg):
            if not isinstance(arg, list) or len(arg) < 2:
                raise ValueError("xdata: must be a 2D list with at least two rows")
            vals = []
            for i, row in enumerate(arg):
                if not isinstance(row, list) or len(row) == 0:
                    raise ValueError(f"xdata row {i}: must be a non-empty list")
                try:
                    vals.append(float(row[0]))
                except Exception:
                    raise ValueError(f"xdata row {i}: non-numeric value")
            return np.asarray(vals, dtype=np.float64)

        def _convert_y(arg):
            if not isinstance(arg, list) or len(arg) < 2:
                raise ValueError("ydata: must be a 2D list with at least two rows")
            rows = []
            expected_len = None
            for i, row in enumerate(arg):
                if not isinstance(row, list) or len(row) == 0:
                    raise ValueError(f"ydata row {i}: must be a non-empty list")
                try:
                    vals = [float(v) for v in row]
                except Exception:
                    raise ValueError(f"ydata row {i}: non-numeric value")
                if expected_len is None:
                    expected_len = len(vals)
                elif len(vals) != expected_len:
                    raise ValueError("ydata rows must all have the same length")
                rows.append(vals)

            y_arr = np.asarray(rows, dtype=np.float64)
            if y_arr.ndim == 1:
                y_arr = y_arr.reshape(-1, 1)
            return y_arr

        x_arr = _convert_x(xdata)
        y_arr = _convert_y(ydata)

        if x_arr.shape[0] != y_arr.shape[0]:
            raise ValueError("xdata and ydata must have the same number of rows")
        return x_arr, y_arr

    # Model definitions dictionary
    models = {
        'two_segment_linear_breakpoint': {
            'params': ['a1', 'k1', 'xi', 'k2'],
            'model': lambda x, a1, k1, xi, k2: np.where(x < xi, a1 + k1 * x, a1 + k1 * xi + k2 * (x - xi)),
            'guess': lambda xa, ya: (float(np.min(ya)), float(np.polyfit(xa, ya, 1)[0]) if xa.size > 1 else 0.0, float(np.median(xa)), float(np.polyfit(xa, ya, 1)[0]) if xa.size > 1 else 0.0),
        },
        'three_segment_linear_breakpoint': {
            'params': ['a1', 'k1', 'xi1', 'k2', 'xi2', 'k3'],
            'model': lambda x, a1, k1, xi1, k2, xi2, k3: np.where(x < xi1, a1 + k1 * x, np.where(x < xi2, a1 + k1 * xi1 + k2 * (x - xi1), a1 + k1 * xi1 + k2 * (xi2 - xi1) + k3 * (x - xi2))),
            'guess': lambda xa, ya: (float(np.min(ya)), 0.0, float(np.percentile(xa, 33)), 0.0, float(np.percentile(xa, 66)), 0.0),
        },
        'heaviside_step_two_level': {
            'params': ['A', 'B', 'x1'],
            'model': lambda x, A, B, x1: np.where(x < x1, A, B),
            'guess': lambda xa, ya: (float(np.percentile(ya, 75)), float(np.percentile(ya, 25)), float(np.median(xa))),
        },
        'langevin_paramagnetic': {
            'params': ['y0', 'xc', 'C'],
            'model': lambda x, y0, xc, C: y0 + C * (np.where(np.abs(x - xc) > 1e-9, np.cosh(x - xc) / np.sinh(x - xc) - 1.0 / (x - xc), 0.0)),
            'guess': lambda xa, ya: (float(np.min(ya)), float(np.median(xa)), float(np.ptp(ya) if np.ptp(ya) else 1.0)),
        },
        'langevin_scaled_field': {
            'params': ['y0', 'xc', 'C', 's'],
            'model': lambda x, y0, xc, C, s: y0 + C * (np.where(np.abs(x - xc) > 1e-9, np.cosh((x - xc) / s) / np.sinh((x - xc) / s) - s / (x - xc), 0.0)),
            'guess': lambda xa, ya: (float(np.min(ya)), float(np.median(xa)), float(np.ptp(ya) if np.ptp(ya) else 1.0), 1.0),
        },
        'planar_surface_two_outputs': {
            'params': ['a', 'b', 'c', 'd'],
            'model': lambda x, a, b, c, d: np.column_stack((a + b * x, c + d * x)).reshape(-1, order='F'),
            'guess': lambda xa, ya: (
                float(np.mean(ya[:, 0])),
                0.0,
                float(np.mean(ya[:, 1])),
                0.0,
            ),
            'output_mode': 'flatten_columns',
        },
        'linear_exponential_multivariate': {
            'params': ['a1', 'a2', 'b1', 'b2', 'c'],
            'model': lambda x, a1, a2, b1, b2, c: np.column_stack((a1 + b1 * x, a2 + b2 * np.exp(c * x))).reshape(-1, order='F'),
            'guess': lambda xa, ya: (
                float(np.mean(ya[:, 0])),
                float(np.mean(ya[:, 1])),
                0.0,
                0.1,
                0.1,
            ),
            'output_mode': 'flatten_columns',
        }
    }

    # Validate model parameter
    if piecewise_model not in models:
        return f"Invalid model: {str(piecewise_model)}. Valid models are: {', '.join(models.keys())}"

    model_info = models[piecewise_model]

    # Validate and convert input data
    try:
        x_arr, y_matrix = _validate_data(xdata, ydata)
    except ValueError as e:
        return f"Invalid input: {e}"

    output_mode = model_info.get('output_mode', 'single')
    if output_mode == 'single':
        if y_matrix.shape[1] != 1:
            return "Invalid input: ydata must contain a single column for the selected model"
        y_arr = y_matrix[:, 0]
        guess_data = y_arr
    elif output_mode == 'flatten_columns':
        if y_matrix.shape[1] < 2:
            return "Invalid input: ydata must contain at least two columns for the selected model"
        y_arr = y_matrix.reshape(-1, order='F')
        guess_data = y_matrix
    else:
        return f"Invalid model configuration: unsupported output mode '{output_mode}'"

    # Perform curve fitting
    try:
        p0 = model_info['guess'](x_arr, guess_data)
        bounds = model_info.get('bounds', (-np.inf, np.inf))
        if bounds == (-np.inf, np.inf):
            popt, pcov = scipy_curve_fit(model_info['model'], x_arr, y_arr, p0=p0, maxfev=10000)
        else:
            popt, pcov = scipy_curve_fit(model_info['model'], x_arr, y_arr, p0=p0, bounds=bounds, maxfev=10000)

        fitted_vals = [float(v) for v in popt]
        for v in fitted_vals:
            if math.isnan(v) or math.isinf(v):
                return "Fitting produced invalid numeric values (NaN or inf)."
    except ValueError as e:
        return f"Initial guess error: {e}"
    except Exception as e:
        return f"curve_fit error: {e}"

    # Calculate standard errors
    std_errors = None
    try:
        if pcov is not None and np.isfinite(pcov).all():
            std_errors = [float(v) for v in np.sqrt(np.diag(pcov))]
    except Exception:
        pass

    return [model_info['params'], fitted_vals, std_errors] if std_errors else [model_info['params'], fitted_vals]

Online Calculator