LOWESS Project

One LOWESS to Rule Them All

The fastest, most robust, and most feature-complete language-agnostic LOWESS (Locally Weighted Scatterplot Smoothing) implementation for Rust, Python, R, Julia, JavaScript, C++, and WebAssembly.

[!IMPORTANT]

The lowess-project contains a complete ecosystem for LOWESS smoothing:

• lowess - Core single-threaded Rust implementation with no_std support
• fastLowess - Parallel CPU and GPU-accelerated Rust wrapper with ndarray integration
• R bindings - extendr-based R binding
• Python bindings - PyO3-based Python binding
• Julia bindings - Native Julia binding with C FFI
• JavaScript bindings - Node.js binding
• WebAssembly bindings - WASM binding
• C++ bindings - Native C++ binding with CMake integration

Installation

[!NOTE]

Currently available for R, Python, Rust, Julia, Node.js, WebAssembly, and C++. See the Installation Guide for detailed installation instructions.

GPU Backend

In addition to parallel = true (multi-core CPU), the batch Lowess class in every binding — except WebAssembly — as well as the fastLowess Rust crate itself, can run on the GPU via wgpu (Vulkan/Metal/DX12). It's opt-in and worth enabling for high-throughput processing of large datasets (roughly 10k+ points); for smaller inputs the CPU backend is typically faster. StreamingLowess/OnlineLowess remain CPU-only. See the GPU Backend guide for installation instructions and usage.

Documentation

[!NOTE]

📚 View the full documentation

LOESS vs. LOWESS

[!TIP] Note: For a LOESS implementation, use loess-project.

Why this package?

Speed

The lowess project beats the competition in terms of speed, whether in single-threaded or multi-threaded parallel execution. It is on average 200-327x faster than Python's statsmodels.lowess and 2-3x faster than R's lowess.

For more details on the performance comparison, see the Benchmarks page.

Robustness

This implementation is more robust than R's lowess and Python's statsmodels due to two key design choices:

MAD-Based Scale Estimation:

For robustness weight calculations, this crate uses Median Absolute Deviation (MAD) for scale estimation:

s = median(|r_i - median(r)|)

In contrast, statsmodels and R's lowess uses the median of absolute residuals (MAR):

s = median(|r_i|)

• MAD is a breakdown-point-optimal estimator—it remains valid even when up to 50% of data are outliers.
• The median-centering step removes asymmetric bias from residual distributions.
• MAD provides consistent outlier detection regardless of whether residuals are centered around zero.

Boundary Padding:

This crate applies a range of different boundary policies at dataset edges:

• Extend: Repeats edge values to maintain local neighborhood size.
• Reflect: Mirrors data symmetrically around boundaries.
• Zero: Pads with zeros (useful for signal processing).
• NoBoundary: Original Cleveland behavior

statsmodels and R's lowess do not apply boundary padding, which can lead to:

• Biased estimates near boundaries due to asymmetric local neighborhoods.
• Increased variance at the edges of the smoothed curve.

Features

A variety of features, supporting a range of use cases:

Validation

All implementations are numerical twins of R's lowess:

API Reference

R:

library(rfastlowess)

model <- Lowess(
    fraction = 0.5,
    iterations = 3L,
    delta = 0.01,
    weight_function = "tricube",
    robustness_method = "bisquare",
    scaling_method = "mad",
    zero_weight_fallback = "use_local_mean",
    boundary_policy = "extend",
    confidence_intervals = 0.95,
    prediction_intervals = 0.95,
    return_diagnostics = TRUE,
    return_residuals = TRUE,
    return_robustness_weights = TRUE,
    return_se = TRUE,
    cv_fractions = c(0.3, 0.5, 0.7),
    cv_method = "kfold",
    cv_k = 5L,
    cv_seed = 123L,
    auto_converge = 1e-4,
    parallel = TRUE
)
custom_weights <- rep(1, length(x))
result <- fit(model, x, y, custom_weights = custom_weights)

# Result structure:
result$x,
result$y,
result$standard_errors,
result$confidence_lower,
result$confidence_upper,
result$prediction_lower,
result$prediction_upper,
result$residuals,
result$robustness_weights,
result$diagnostics,
result$iterations_used,
result$fraction_used,
result$cv_scores

Python:

from fastlowess import Lowess

model = Lowess(
    fraction=0.5,
    iterations=3,
    delta=0.01,
    weight_function="tricube",
    robustness_method="bisquare",
    scaling_method="mad",
    zero_weight_fallback="use_local_mean",
    boundary_policy="extend",
    confidence_intervals=0.95,
    prediction_intervals=0.95,
    return_diagnostics=True,
    return_residuals=True,
    return_robustness_weights=True,
    return_se=True,
    cv_fractions=[0.3, 0.5, 0.7],
    cv_method="kfold",
    cv_k=5,
    cv_seed=123,
    auto_converge=1e-4,
    parallel=True
)
custom_weights = [1.0] * len(x)
result = model.fit(x, y, custom_weights=custom_weights)

# Result structure:
result.x,
result.y,
result.standard_errors,
result.confidence_lower,
result.confidence_upper,
result.prediction_lower,
result.prediction_upper,
result.residuals,
result.robustness_weights,
result.diagnostics,
result.iterations_used,
result.fraction_used,
result.cv_scores

Rust:

use lowess::prelude::*;

let model = Lowess::new()
    .fraction(0.5)
    .iterations(3)
    .delta(0.01)
    .weight_function("tricube")
    .robustness_method("bisquare")
    .scaling_method("mad")
    .zero_weight_fallback("use_local_mean")
    .boundary_policy("extend")
    .return_se()
    .confidence_intervals(0.95)
    .prediction_intervals(0.95)
    .return_diagnostics()
    .return_residuals()
    .return_robustness_weights()
    .cv_method("kfold")
    .cv_k(5)
    .cv_fractions(vec![0.3, 0.5, 0.7])
    .cv_seed(123)
    .auto_converge(1e-4)
    .custom_weights(vec![1.0; x.len()])
    .build()?;

let result = model.fit(&x, &y)?;

// Result structure:
pub struct LowessResult<T> {
    pub x: Vec<T>,                           // Sorted x values
    pub y: Vec<T>,                           // Smoothed y values
    pub standard_errors: Option<Vec<T>>,
    pub confidence_lower: Option<Vec<T>>,
    pub confidence_upper: Option<Vec<T>>,
    pub prediction_lower: Option<Vec<T>>,
    pub prediction_upper: Option<Vec<T>>,
    pub residuals: Option<Vec<T>>,
    pub robustness_weights: Option<Vec<T>>,
    pub diagnostics: Option<Diagnostics<T>>,
    pub iterations_used: Option<usize>,
    pub fraction_used: T,
    pub cv_scores: Option<Vec<T>>,
}

Julia:

using FastLOWESS

model = Lowess(;
    fraction=0.5,
    iterations=3,
    delta=NaN,  # NaN for auto
    weight_function="tricube",
    robustness_method="bisquare",
    scaling_method="mad",
    zero_weight_fallback="use_local_mean",
    boundary_policy="extend",
    confidence_intervals=NaN,
    prediction_intervals=NaN,
    return_diagnostics=true,
    return_residuals=true,
    return_robustness_weights=true,
    return_se=true,
    cv_fractions=Float64[], # e.g. [0.3, 0.5]
    cv_method="kfold",
    cv_k=5,
    cv_seed=123,
    auto_converge=NaN,
    parallel=true
)
custom_weights = ones(length(x))
result = fit(model, x, y; custom_weights=custom_weights)

# Result structure:
result.x,
result.y,
result.standard_errors,
result.confidence_lower,
result.confidence_upper,
result.prediction_lower,
result.prediction_upper,
result.residuals,
result.robustness_weights,
result.diagnostics,
result.iterations_used,
result.fraction_used,
result.cv_scores

Node.js:

import { Lowess } from "fastlowess"

const model = new Lowess({
    fraction: 0.5,
    iterations: 3,
    delta: 0.01,
    weight_function: "tricube",
    robustness_method: "bisquare",
    scaling_method: "mad",
    zero_weight_fallback: "use_local_mean",
    boundary_policy: "extend",
    return_se: true,
    confidence_intervals: 0.95,
    prediction_intervals: 0.95,
    return_diagnostics: true,
    return_residuals: true,
    return_robustness_weights: true,
    cv_fractions: [0.3, 0.5, 0.7],
    cv_method: "kfold",
    cv_k: 5,
    cv_seed: 123,
    auto_converge: 1e-4,
    parallel: true
})
const custom_weights = Array(x.length).fill(1.0)
const result = model.fit(x, y, custom_weights)

// Result structure:
result.x,
result.y,
result.standard_errors,
result.confidence_lower,
result.confidence_upper,
result.prediction_lower,
result.prediction_upper,
result.residuals,
result.robustness_weights,
result.diagnostics,
result.iterations_used,
result.fraction_used,
result.cv_scores

WebAssembly:

import { Lowess } from "fastlowess-wasm"

const model = new Lowess({
    fraction: 0.5,
    iterations: 3,
    delta: 0.01,
    weight_function: "tricube",
    robustness_method: "bisquare",
    scaling_method: "mad",
    zero_weight_fallback: "use_local_mean",
    boundary_policy: "extend",
    return_se: true,
    confidence_intervals: 0.95,
    prediction_intervals: 0.95,
    return_diagnostics: true,
    return_residuals: true,
    return_robustness_weights: true,
    cv_fractions: [0.3, 0.5, 0.7],
    cv_method: "kfold",
    cv_k: 5,
    cv_seed: 123,
    auto_converge: 1e-4,
    parallel: true
})
const custom_weights = new Float64Array(x.length).fill(1)
const result = model.fit(x, y, custom_weights)

// Result structure:
result.x,
result.y,
result.standard_errors,
result.confidence_lower,
result.confidence_upper,
result.prediction_lower,
result.prediction_upper,
result.residuals,
result.robustness_weights,
result.diagnostics,
result.iterations_used,
result.fraction_used,
result.cv_scores

C++:

#include "fastlowess.hpp"

fastlowess::LowessOptions options;
options.fraction = 0.5;
options.iterations = 3;
options.delta = 0.01;
options.weight_function = "tricube";
options.robustness_method = "bisquare";
options.scaling_method = "mad";
options.zero_weight_fallback = "use_local_mean";
options.boundary_policy = "extend";
options.confidence_intervals = 0.95;
options.prediction_intervals = 0.95;
options.return_diagnostics = true;
options.return_residuals = true;
options.return_robustness_weights = true;
options.return_se = true;
options.cv_fractions = {0.3, 0.5, 0.7};
options.cv_method = "kfold";
options.cv_k = 5;
options.cv_seed = 123;
options.auto_converge = 1e-4;
options.parallel = true;

fastlowess::Lowess model(options);
std::vector<double> custom_weights(x.size(), 1.0);
const auto result = model.fit(x, y, custom_weights).value();

// Result structure:
result.x_vector(),
result.y_vector(),
result.standard_errors(),
result.confidence_lower(),
result.confidence_upper(),
result.prediction_lower(),
result.prediction_upper(),
result.residuals(),
result.robustness_weights(),
result.diagnostics(),
result.iterations_used(),
result.fraction_used(),
result.cv_scores()

Contributing

Contributions are welcome! Please see the Contributing Guide for more information.

Changelog

See the Changelog for a history of changes.

License

Licensed under MIT or Apache-2.0.

Citation

If you use this software in your research, please cite it using the CITATION.cff file or the BibTeX entry below:

@software{lowess_project,
  author = {Valizadeh, Amir},
  title = {LOWESS Project: High-Performance Locally Weighted Scatterplot Smoothing},
  year = {2026},
  url = {https://github.com/thisisamirv/lowess-project},
  license = {MIT OR Apache-2.0}
}