Add JOSS 2.0 paper draft

.github/workflows/draft-pdf.yml

@@ -0,0 +1,19 @@
+on: [push]
+
+jobs:
+ paper:
+ runs-on: ubuntu-latest
+ name: Paper Draft
+ steps:
+ - name: Checkout
+ uses: actions/checkout@v6
+ - name: Build draft PDF
+ uses: openjournals/openjournals-draft-action@master
+ with:
+ journal: joss
+ paper-path: paper/paper.md
+ - name: Upload
+ uses: actions/upload-artifact@v7
+ with:
+ name: paper
+ path: paper/paper.pdf

LICENSE.txt

@@ -1,5 +1,3 @@
-The MIT License (MIT)
-
 Copyright (c) 2025
 Floris van Breugel, The University of Nevada, Reno, Mechanical Engineering Department
 Pavel Komarov, Univerity of Washington Electrical and Computer Engineering

paper/paper.bib

@@ -0,0 +1,167 @@
+@article{vanBreugel2022,
+ doi = {10.21105/joss.04078},
+ url = {https://doi.org/10.21105/joss.04078},
+ year = {2022},
+ publisher = {The Open Journal},
+ volume = {7},
+ number = {71},
+ pages = {4078},
+ author = {{van Breugel}, Floris and Liu, Yuying and Brunton, Bingni W. and Kutz, J. Nathan},
+ title = {{PyNumDiff}: A {Python} package for numerical differentiation of noisy time-series data},
+ journal = {Journal of Open Source Software}
+}
+
+@misc{komarov2025,
+ title = {A Taxonomy of Numerical Differentiation Methods},
+ author = {Komarov, Pavel and {van Breugel}, Floris and Kutz, J. Nathan},
+ year = {2025},
+ eprint = {2512.09090},
+ archivePrefix = {arXiv},
+ primaryClass = {math.NA},
+ url = {https://arxiv.org/abs/2512.09090}
+}
+
+@article{vanBreugel2020numerical,
+ doi = {10.1109/ACCESS.2020.3034077},
+ year = {2020},
+ author = {{van Breugel}, Floris and Kutz, J. Nathan and Brunton, Bingni W.},
+ journal = {IEEE Access},
+ title = {Numerical differentiation of noisy data: A unifying multi-objective optimization framework},
+ volume = {8},
+ pages = {196865--196877}
+}
+
+@article{chartrand2011numerical,
+ author = {Rick Chartrand},
+ title = {Numerical differentiation of noisy, nonsmooth data},
+ journal = {ISRN Applied Mathematics},
+ year = {2011},
+ volume = {2011},
+ pages = {164564},
+ doi = {10.5402/2011/164564}
+}
+
+@article{brunton2016discovering,
+ author = {Steven L. Brunton and Joshua L. Proctor and J. Nathan Kutz},
+ title = {Discovering governing equations from data by sparse identification of nonlinear dynamical systems},
+ journal = {Proceedings of the National Academy of Sciences},
+ year = {2016},
+ volume = {113},
+ number = {15},
+ pages = {3932--3937},
+ doi = {10.1073/pnas.1517384113}
+}
+
+@article{virtanen2020scipy,
+ author = {Pauli Virtanen and Ralf Gommers and Travis E. Oliphant and others},
+ title = {{SciPy} 1.0: Fundamental algorithms for scientific computing in {Python}},
+ journal = {Nature Methods},
+ year = {2020},
+ volume = {17},
+ pages = {261--272},
+ doi = {10.1038/s41592-019-0686-2}
+}
+
+@article{harris2020array,
+ author = {Charles R. Harris and K. Jarrod Millman and St{\'{e}}fan J. {van der Walt} and others},
+ title = {Array programming with {NumPy}},
+ journal = {Nature},
+ year = {2020},
+ volume = {585},
+ pages = {357--362},
+ doi = {10.1038/s41586-020-2649-2}
+}
+
+@article{diamond2016cvxpy,
+ author = {Steven Diamond and Stephen Boyd},
+ title = {{CVXPY}: A {Python}-embedded modeling language for convex optimization},
+ journal = {Journal of Machine Learning Research},
+ year = {2016},
+ volume = {17},
+ number = {83},
+ pages = {1--5}
+}
+
+@article{savitzky1964,
+ author = {Abraham Savitzky and Marcel J. E. Golay},
+ title = {Smoothing and Differentiation of Data by Simplified Least Squares Procedures},
+ journal = {Analytical Chemistry},
+ year = {1964},
+ volume = {36},
+ number = {8},
+ pages = {1627--1639},
+ doi = {10.1021/ac60214a047}
+}
+
+@article{rauch1965,
+ author = {Herbert E. Rauch and F. Tung and Charlotte T. Striebel},
+ title = {Maximum likelihood estimates of linear dynamic systems},
+ journal = {AIAA Journal},
+ year = {1965},
+ volume = {3},
+ number = {8},
+ pages = {1445--1450},
+ doi = {10.2514/3.3166}
+}
+
+@article{aravkin2013,
+ author = {Aleksandr Y. Aravkin and James V. Burke and Gianluigi Pillonetto},
+ title = {Optimization viewpoint on {Kalman} smoothing with applications to robust and sparse estimation},
+ journal = {Journal of Machine Learning Research},
+ year = {2013},
+ volume = {14},
+ pages = {2513--2558},
+ url = {https://jmlr.org/papers/volume14/aravkin13a/aravkin13a.pdf},
+ doi = {10.1007/978-3-642-38398-4_8}
+}
+
+@software{derivative_pkg,
+ author = {Andy Goldschmidt},
+ title = {derivative: Numerical differentiation in {Python}},
+ year = {2021},
+ url = {https://github.com/andgoldschmidt/derivative}
+}
+
+@article{kalman1960,
+ author = {Rudolf E. Kalman},
+ title = {A new approach to linear filtering and prediction problems},
+ journal = {Journal of Basic Engineering},
+ year = {1960},
+ volume = {82},
+ number = {1},
+ pages = {35--45},
+ doi = {10.1115/1.3662552}
+}
+
+@article{ahnert2007,
+ author = {Karsten Ahnert and Markus Abel},
+ title = {Numerical differentiation of experimental data: local versus global methods},
+ journal = {Computer Physics Communications},
+ year = {2007},
+ volume = {177},
+ number = {10},
+ pages = {764--774},
+ doi = {10.1016/j.cpc.2007.03.009}
+}
+
+@article{pysindy,
+ author = {Brian M. de Silva and Kathleen Champion and Markus Quade and Jean-Christophe Loiseau and J. Nathan Kutz and Steven L. Brunton},
+ title = {{PySINDy}: A {Python} package for the sparse identification of nonlinear dynamics from data},
+ journal = {Journal of Open Source Software},
+ year = {2020},
+ volume = {5},
+ number = {49},
+ pages = {2104},
+ doi = {10.21105/joss.02104}
+}
+
+@article{lee2019pywavelets,
+ author = {Gregory R. Lee and Ralf Gommers and Filip Wasilewski and Kai Wohlfahrt and Aaron O'Leary},
+ title = {{PyWavelets}: A {Python} package for wavelet analysis},
+ journal = {Journal of Open Source Software},
+ year = {2019},
+ volume = {4},
+ number = {36},
+ pages = {1237},
+ doi = {10.21105/joss.01237}
+}

paper/paper.md

@@ -0,0 +1,109 @@
+---
+title: 'PyNumDiff: Practical Numerical Differentiation for Noisy Data'
+tags:
+ - Python
+ - numerical differentiation
+ - time series
+ - denoising
+ - dynamics
+ - signal processing
+authors:
+ - name: Pavel Komarov
+ orcid: 0009-0007-7482-2807
+ corresponding: true
+ affiliation: 1
+ - name: Floris van Breugel
+ orcid: 0000-0001-6538-7179
+ affiliation: 2
+ - name: Maria Protogerou
+ affiliation: 4
+ - name: J. Nathan Kutz
+ orcid: 0000-0002-6004-2275
+ affiliation: 3
+affiliations:
+ - name: Department of Electrical and Computer Engineering, University of Washington, USA
+ index: 1
+ - name: Department of Mechanical Engineering, University of Nevada, Reno, USA
+ index: 2
+ - name: Autodesk Research, London, UK
+ index: 3
+ - name: Department of Applied Mathematics, University of Washington, USA
+ index: 4
+date: 1 April 2026
+bibliography: paper.bib
+---
+
+# Summary
+
+Derivatives of measured data are a prerequisite across science and engineering: identifying governing equations, designing controllers, and processing sensor streams alike. The textbook remedy, finite differencing, amplifies noise as $1/\Delta t$ and deteriorates rapidly as data grows noisier or more finely sampled. Smoothing before differencing helps, but algorithm choice and tuning substantially affect the result, and no single approach wins universally.
+
+PyNumDiff is an open-source Python package consolidating a broad suite of numerical differentiation methods under a unified API. Seven algorithm families are implemented: (1) prefiltering followed by finite difference calculation; (2) iterated finite differencing; (3) polynomial fitting [@savitzky1964]; (4) spectral, radial basis function, and wavelet fitting; (5) total variation regularization [@chartrand2011numerical]; (6) Kalman smoothing [@kalman1960; @rauch1965]; and (7) local approximation with linear models. @ahnert2007 provides a useful taxonomy of these families, distinguishing local methods (which estimate derivatives from a surrounding window) from global methods that fit the entire signal at once. All PyNumDiff methods return a matched pair `(x_hat, dxdt_hat)`. A companion paper [@komarov2025] benchmarks all methods across test signals and guides selection for different application scenarios. This paper describes the package's second generation, the 0.2.x release series.
+
+
+# Statement of Need
+
+Numerical differentiation has a diverse ecosystem of specialized algorithms, each with different strengths in outlier robustness, computational cost, irregular sampling, or missing observations, but no consolidated home. Without one, practitioners are left assembling solutions piecemeal from disparate packages [@vanBreugel2022].
+
+PyNumDiff addresses this gap. Its unified interface lets users compare methods on the same data, exploit specialized capabilities, and select hyperparameters without ground-truth derivatives. Derivative estimation trades data fidelity against smoothness, and PyNumDiff frames that tradeoff as a multi-objective optimization spanning all methods, with a single scalar weight steering between the two. When a ground-truth derivative is available, it can serve as the optimization target instead [@vanBreugel2020numerical]. A natural and growing application is SINDy [@brunton2016discovering], which discovers governing equations by regressing measured derivatives, making the quality of those estimates a direct determinant of model accuracy.
+
+
+# State of the Field
+
+Relevant Python tools exist, but none covers PyNumDiff's breadth. `numpy.gradient` and `scipy.signal.savgol_filter` [@virtanen2020scipy] handle only a sliver of the method space; `findiff` offers high-order finite difference stencils suited to clean simulation data, not noisy measurements. Historically, practitioners have had to stitch together PyKalman, PyDMD, and standalone TVR scripts [@chartrand2011numerical] with no shared API or principled way to compare results. The `derivative` package [@derivative_pkg] overlaps substantially but lacks multidimensional support, NaN handling, and hyperparameter optimization. No existing package spans PyNumDiff's seven method families with a consistent interface.
+
+The original PyNumDiff publication [@vanBreugel2022] established the core method set and optimization framework. This second generation rewrites it from the ground up and consolidates it. `kerneldiff`, `rtsdiff`, and `tvrdiff` each cover by parameter what had been separate functions, so the package presents fewer entry points than before while spanning more methods. `rbfdiff`, `waveletdiff`, and `robustdiff` are new, bringing the set to the twelve methods of Table 1, all organized into the seven families above behind one keyword-argument signature. A single optimizer serves every one of them, and four capabilities now span the package, covering multidimensional data, irregular sample spacing, missing observations, and circular domains.
+
+
+# Software Design
+
+**Package design.** All differentiation methods share the call signature
+
+```python
+x_hat, dxdt_hat = method(x, dt_or_t, **params)
+```
+
+where `x` is a NumPy array [@harris2020array] of measurements; `dt_or_t` is either a scalar step size or an array of sample locations; and keyword arguments configure the method, making calls self-documenting. Prior positional signatures are preserved with deprecation warnings.
+
+**Software architecture.** PyNumDiff is organized into seven method modules plus shared `utils` and `optimize` modules in a flat structure. Where strong alternatives exist, PyNumDiff delegates rather than reimplements: SciPy [@virtanen2020scipy] provides spline fitting, Savitzky-Golay filtering, and signal processing routines; NumPy [@harris2020array] provides the FFT; PyWavelets [@lee2019pywavelets] provides the discrete wavelet transform for `waveletdiff`; CVXPY [@diamond2016cvxpy] handles convex optimization for `robustdiff` and `tvrdiff` as an optional dependency. The public `kalman_filter` and `rts_smooth` primitives let users with known dynamics bypass `rtsdiff`'s constant-derivative model.
+
+**Method capabilities.** Table 1 groups the twelve non-deprecated methods by capability. All handle multidimensional data via `axis`, and several go further.
+
+| Capability | Methods |
+|---|---|
+| Multidimensional data | `kerneldiff`, `finitediff`, `polydiff`, `savgoldiff`, `splinediff`, `spectraldiff`, `rbfdiff`, `waveletdiff`, `tvrdiff`, `rtsdiff`, `robustdiff`, `lineardiff` |
+| Variable step | `polydiff`, `splinediff`, `rbfdiff`, `rtsdiff`, `robustdiff` |
+| Missing data | `polydiff`, `splinediff`, `rtsdiff`, `robustdiff` |
+| Outlier robustness | `tvrdiff`, `robustdiff` |
+| Circular domain | `rtsdiff` |
+
+Table: Methods by capability.
+
+**Irregular and incomplete sampling.** Methods that support variable step size accept an array of sample locations in place of a scalar step, and Kalman-based methods then compute the transition matrix by matrix exponential at each actual interval. NaN entries are treated as missing observations, excluded from fitting and imputed from the model, so sensors that drop samples need no preprocessing.
+
+**Outlier robustness.** `robustdiff` replaces the quadratic Kalman cost with Huber loss terms on both measurement and process residuals, following @aravkin2013, with CVXPY [@diamond2016cvxpy] as the optimization backend; the sparse problem formulation scales linearly with signal length. `tvrdiff` similarly applies Huber loss on data fidelity; its total variation penalty on the derivative additionally promotes piecewise-smooth solutions for signals with abrupt transitions.
+
+**Circular and wrapped domains.** `rtsdiff` accepts `circular=True` for quantities like angles on a periodic domain. Innovation residuals are wrapped to $[-\pi, \pi]$ before each Kalman update via an `innovation_fn` hook, and `x_hat` is returned in the same range, avoiding the large spurious spikes naive smoothers produce when a signal crosses the $\pm\pi$ boundary.
+
+**Hyperparameter optimization.** `pynumdiff.optimize` selects hyperparameters $\Phi$ by minimizing
+$$L(\Phi) = \text{RMSE}\big(\textstyle\int\hat{\dot{x}}(\Phi) + c,\; x\big) + \gamma\,\text{TV}\big(\hat{\dot{x}}(\Phi)\big),$$
+which requires no ground truth, because fidelity is measured by integrating the estimated derivative back against the measured signal [@vanBreugel2020numerical]. The smoothness weight $\gamma$, exposed as `tvgamma`, can be initialized automatically from the signal's estimated cutoff frequency. This version robustifies the first term with a Huber penalty so outliers do not bias parameter selection, and it reduces the Kalman parameter space from two independent noise variances to their log-ratio, the only salient factor [@komarov2025]. Categorical and boolean hyperparameters are now supported natively, so discrete choices such as derivative order can be optimized jointly with continuous ones. Repeated evaluations are cached to avoid redundant work. A single call to `suggest_method` runs this search across nearly every method in the package and reports the best fit, which is the practical payoff of putting them all behind one interface.
+
+**Testing and continuous integration.** The test suite validates all methods against analytic functions with known derivatives, covering both noiseless and noisy cases. Care was taken to avoid tautological tests where the implementation directly determines the expected result. GitHub Actions runs the suite on every push and pull request, and Coveralls tracks line coverage, currently 90%.
+
+
+# Research Impact Statement
+
+The original PyNumDiff paper [@vanBreugel2022] has been applied in experimental biology (flight kinematics from motion capture), control engineering (observer design), and data-driven dynamics identification via SINDy [@brunton2016discovering]. Tutorial notebooks in the repository cover basic usage, hyperparameter optimization, automatic method suggestion, a head-to-head performance analysis across simulations, and demonstrations of the outlier, multidimensional, and circular-domain capabilities described above. Each runs end to end against the released package, so the claims made here can be reproduced directly. Full API documentation is published at [pynumdiff.readthedocs.io](https://pynumdiff.readthedocs.io/master/). The companion Taxonomy paper [@komarov2025], submitted to the Journal of Computational Physics, provides the theoretical underpinning and benchmarks all included methods. The PySINDy project [@pysindy] maintains its own differentiation submodule substantially overlapping with PyNumDiff's capabilities; integration discussions are ongoing.
+
+
+# AI Usage Disclosure
+
+This paper was drafted with assistance from Claude Sonnet 4.6 and Claude Opus 4.8 (Anthropic), which also implemented successive code revisions to address recent issues and author feedback. All outputs were reviewed and further edited by hand, and the authors take full responsibility for accuracy.
+
+
+# Acknowledgements
+
+The authors thank Yuying Liu and Bingni W. Brunton for their contributions to the original PyNumDiff package [@vanBreugel2022], and Sasha Aravkin for discussions on convex optimization techniques that informed the robust differentiation methods. This work was supported by the NSF AI Institute in Dynamic Systems (grant number 2112085).
+
+
+# References