Skip to content

Latest commit

 

History

6 Commits

Folders and files

NameName
Last commit message
Last commit date
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

Repository files navigation

LMomFit: L-moment distribution fitting for small samples containing extremes

DOI PyPI View on File Exchange

LMomFit fits probability distributions to small samples that may contain genuine extreme observations, in parallel MATLAB and Python implementations.

Small samples with an extreme value are the awkward case. Conventional moment estimators involve third and fourth powers of deviations, so a single extreme dominates them — and the usual remedy, declaring the observation an outlier and discarding it, throws away exactly the tail information the analysis exists to characterise. LMomFit uses L-moments instead: expectations of linear combinations of order statistics, which exist whenever the mean exists, are nearly unbiased at small n, and are far less sensitive to extremes.

The pipeline is four steps, available as a single guarded call:

  1. Identify a candidate family from the L-moment ratio diagram (L-skewness vs. L-kurtosis).
  2. Estimate its parameters in closed form — no iterative optimisation, so the fit is cheap enough to run inside an outer propagation or design loop.
  3. Evaluate the fitted PDF/CDF, draw variates, and score the fit against a conventional-moment (maximum-likelihood) alternative.
  4. Report honestly — bootstrap selection frequencies with an ambiguity flag, rather than a single family that a scarce sample cannot justify.

Two components have no counterpart in existing L-moment software: a guarded ranked-fallback estimator that makes unattended automatic fitting safe inside batch jobs and optimisation loops, and a bootstrap identification routine that reports a candidate set instead of a single verdict.

Background and validation on statistical distributions and engineering case studies (sheet-metal forming, speed reducer design, probabilistic fatigue life) are described in:

Jayaraman D, Ramu P. L-moments-based uncertainty quantification for scarce samples including extremes. Structural and Multidisciplinary Optimization. 2021 Aug;64(2):505–39.


Installation

Python

pip install lmomfit            # core: NumPy + SciPy only
pip install "lmomfit[ui]"      # adds Streamlit, matplotlib, pandas

The [ui] extra installs the dependencies used by the Streamlit app and the demo script. Those two are repository scripts, not part of the installed package, so run them from a clone of this repository.

Requires Python ≥ 3.9. Imports as lmomfit — a namespace distinct from the unrelated lmoments and lmoments3 packages, and from lmfit, an unrelated least-squares minimisation package with a similar name.

MATLAB

No installation beyond adding the folder to your path:

addpath('path/to/LMomFit')

No toolboxes required. The toolbox depends only on functions core to both MATLAB and Octave — luq_dist.m supplies closed-form PDF/CDF/inverse-CDF for all ten families, so the Statistics and Machine Learning Toolbox is no longer needed. The one exception is the optional demo_example.m, which compares against fitdist; demo_octave.m is the portable equivalent.

Developed and tested on MATLAB R2018b or later.

GNU Octave

pkg install octave/lmomfit-1.0.0.tar.gz
pkg load lmomfit

Rebuild the tarball after changing any .m source with python octave/build_octave_package.py; it assembles the package from the sources in this directory, so the MATLAB and Octave versions cannot drift apart. Verified on GNU Octave 11.3.0, where tests/octave_verify.m runs with no packages loaded and all checks pass.

On Windows, pkg install needs Octave's own usr\bin on the PATH ahead of C:\Windows\System32, otherwise pkg picks up the system tar.exe, which cannot read the MSYS-style paths Octave passes it.


Quick start

Python

import numpy as np
from lmomfit import fit_best

# A scarce sample containing one genuine extreme
x = np.array([197.8, 224.7, 233.8, 236.8, 243.1, 266.9,
              270.3, 285.2, 291.7, 409.0, 547.0, 758.8])

fit = fit_best(x)          # identify + estimate, guarded, in one call
fit.summary()

fit.distribution           # 'generalized extreme value'
fit.parameters             # closed-form estimate
fit.used_fallback          # True if the nearest family was skipped
fit.skipped                # [(family, reason), ...]

fit.cdf(600)               # evaluate the fitted distribution
fit.ppf(0.99)              # 99th percentile
fit.interval(0.95)         # central 95% interval
fit.rvs(1000)              # draw variates

fit.js_div(x)              # score the fit against data
fit.ks_stat(x)             # binning-free alternative

fit_best returns an LMomentFit; identify_dist, identify_dist_bootstrap and parameter_identify likewise return IdentificationResult, BootstrapIdentification and CandidateFits. All of them keep 1.x mapping access (fit["distribution"]), so existing code and replication scripts run unchanged.

When the sample is scarce, do not trust a single reported family:

from lmomfit import identify_dist_bootstrap

boot = identify_dist_bootstrap(x, n_boot=1000)
boot.status                    # 'clear' or 'ambiguous'
boot.selection_frequencies     # [(family, frequency), ...]
boot.t3_ci, boot.t4_ci         # percentile intervals

MATLAB / Octave

% A scarce sample with one extreme value
X = Random_l('lognormal', [0, 0.5, 0], 12, 1);
X(end+1) = 8*max(X);

% Identify and estimate, skipping any family whose closed-form
% estimator is undefined for this sample
[Distribution, Parameter, skipped] = fit_best(X);

% Evaluate the fit
grid     = linspace(0, max(X), 200);
pdf_vals = PDF_l(grid, Distribution, Parameter);
cdf_vals = CDF_l(grid, Distribution, Parameter);

% Score it against a larger reference sample -- no manual binning needed
X_reference = Random_l('lognormal', [0, 0.5, 0], 20000, 1);
JSDiv(Distribution, Parameter, X_reference)
KSStat(Distribution, Parameter, X_reference)

See demo_octave.m for a complete walkthrough that runs in both MATLAB and Octave, or demo_example.m for the MATLAB-only version that also compares against a fitdist fit and plots both.


Function reference

The two languages share a vocabulary, differing only in capitalisation.

Python MATLAB Purpose
lmom lmom.m First n sample L-moments of a data vector.
pwm, l_moment_ratios Underlying probability-weighted moments; ratios as a mapping.
identify_dist Identify_dist.m Rank families by position on the L-moment ratio diagram.
identify_dist_bootstrap Identify_dist_bootstrap.m Uncertainty-aware identification: selection frequencies across bootstrap resamples, percentile intervals for (t3, t4), and an ambiguity flag.
parameter_estimation Parameter_estimation.m Closed-form parameter estimate from sample L-moments, under explicit domain guards.
parameter_identify parameter_identify.m Fit the top k feasible candidates, so families that sit close together on the ratio diagram can be compared rather than collapsed to the nearest.
fit_best fit_best.m Preferred entry point. Identify and fit in one call, walking past any family whose estimator is undefined and reporting what it skipped.
pdf_l, cdf_l, random_l PDF_l.m, CDF_l.m, Random_l.m Evaluate density/CDF, generate variates.
kl_div, js_div KLDiv.m, JSDiv.m Kullback–Leibler / Jensen–Shannon divergence. Takes either two binned mass vectors, or a fit and a raw sample, binning internally.
ks_stat, ks_test KSStat.m Kolmogorov–Smirnov statistic between a fit and a raw sample. Binning-free, so it can confirm a divergence comparison is not a bin-width artefact.
luq_dist.m Closed-form PDF/CDF/inverse-CDF for all ten families. The reason no Statistics toolbox is needed.
LegendreShiftPoly.m, luq_bin_fit.m, luq_percentile.m Internal helpers.

Supported distribution families

uniform, normal, exponential, gumbel, logistic, generalized extreme value, generalized pareto, lognormal, gamma.

Weibull: supported explicitly, excluded from auto-identification

The three-parameter Weibull (weibul) is fully supported when requested by name — estimation, PDF, CDF and sampling all handle it. It is deliberately not among the families the automatic search selects: its L-moment ratio curve passes through or near other families' loci (shape k=1 is the exponential point; near k≈3.6 it sits essentially on the normal point), so including it makes identification ambiguous rather than better. Both implementations behave identically here.


Correlated Latin Hypercube sampling (external dependency)

Earlier versions of this repository bundled lhsgeneral.m (correlated Latin Hypercube sampling) by Iman Moazzen (2060 Project, IESVic, University of Victoria, BC, Canada). It has been removed because its redistribution licence could not be confirmed, and it is not part of this toolbox's core identify/estimate/evaluate pipeline. If your workflow needs the correlated-sampling step described in the companion papers, obtain it from the original author's MATLAB File Exchange entry: https://www.mathworks.com/matlabcentral/fileexchange/56384-lhsgeneral-pd-correlation-n. Everything remaining in this repository is under the MIT licence.

Citing this software

Please cite both the software and the method paper. Machine-readable metadata is in CITATION.cff.

Jayaraman D, Ramu P. L-moments-based uncertainty quantification for scarce samples including extremes. Structural and Multidisciplinary Optimization. 2021 Aug;64(2):505–39. https://doi.org/10.1007/s00158-021-02930-2

License

MIT — see LICENSE, with a carve-out for the third-party file noted above.

Support

Questions and issues: deepanjayram@gmail.com or via the GitHub issue tracker.

Issues and pull requests are welcome.

About

L-moment distribution fitting for small samples containing extremes. Parallel MATLAB and Python, with guarded automatic fitting and bootstrap family selection.

Topics

Resources

Stars

0 stars

Watchers

0 watching

Forks

Releases

Packages

Contributors

Languages