Software:
I'd like to briefly point out an issue that I have experienced a few times. I
think it's important for scientific practitioners to clean up and package code
to the point where other people who are interested in the work have a fighting
chance of actually running the code themselves. I do my best to meet that
expectation. On the other hand, however, I am not a software engineer, and so my
code is not always perfect in terms of steering users away from footguns or
testing all edge cases. So if you try any of this code, which I sincerely hope
you do, and experience something working way differently (often much worse) than
a paper or figure or example file suggests, please reach out to me before
writing off the code or method completely. There is a non-negligible chance that
I just haven't enforced some aspect of proper usage that may be obvious to a
practitioner in that field---such as sensibly enumerating the observations when
assembling a hierarchical matrix approximation---but that may not occur to
somebody who is just trying out the code or who isn't particularly familiar with
the field. While I obviously do my best to try and enforce proper usage through
the code design itself, it is challenging and I make mistakes. So please do
reach out if you experience surprising things when trying any of these software
packages. There might be an easy remedy, and you might save future users similar
frustrations.
As a second small rant, if you use any of these packages in your own research in
a way that substantively helped your project, *please* do cite them
appropriately. Almost all of these packages are companions to a paper, or the
README of the package will have citation instructions. As an example, please
cite the "Fitting Matern..." paper if you use BesselK.jl to fit a full
three-parameter Matern model. It may seem like sort of a trivial software
dependency, but in this hypothetical you're really using that package to do
something that for decades hasn't really been reliably possible. So please do
acknowledge the work that did make your results possible.
I write all software using the GPLv2 license unless otherwise specified.
-
SMultitaper.jl,
a software package that provides a simple but efficient implementation of
an adaptively weighted multitaper spectral density estimator. There is a
lot of workspace and result caching, so even without careful
pre-allocations repeated calls will be fast. You'll see that this package
doesn't have many commits in the last few years, but that's because it's
done and it has basically no dependencies that need to be bumped. The code
still works well and I do still use it regularly. If you want to compute
hundreds of adaptively weighted multitaper estimates for data of the same
(padded) length, this package is probably a plenty good choice.
-
GPMaxlik.jl,
a simple software package that provides a reasonably tuned function for
evaluating the exact log-likelihood for Gaussian processes and stochastic
estimators for its first and second derivatives. It also now provides
exact methods for the log-likelihood where I have added custom methods for
the gradient and Hessian, so that you can call `ForwardDiff.hessian` on it
and it will actually compute the Hessian using the analytical formula,
thus using LAPACK and only operate on arrays of doubles, but it will still
use ForwardDiff for the derivatives of the kernel function.
-
HalfSpectral.jl,
a sort of software companion to "Flexible nonstationary..." above, this
small package provides a reasonably tuned method for conveniently working
with half-spectral covariance functions for regular data completely in the
time domain. The example files demonstrate how straightforward it is to
build quite complicated models with very little code and then use
automatic differentiation for the derivatives. Besides the recent re-write
for complete AD compatibility, the example files also provide a
demonstration of how to efficiently hook it in to `Vecchia.jl` so that you
can fit your fancy half-spectral covariance function scalably and
efficiently.
-
Vecchia.jl,
(MIT license) A very simple implementation of Vecchia's likelihood
approximation applied to mean-zero Gaussian processes. The code is
specifically optimized to have low allocations and perform well with
multiple threads (in particular, by using aggressive pre-allocation to
squeeze out all allocations from the section that spawns tasks to compute
the small likelihood terms). Gradients and true Hessians are available
using automatic differentiation (ForwardDiff and ReverseDiff), and are
very fast. The package also offers the "reverse" Cholesky factor (a la
Katzfuss and Guinness 2021), and also gives a convenient software
implementation of the EM-based method for estimating parameters when you
have measurement noise that was described in "A Scalable Method to Exploit
Screening...". See the README for a few more tricks that it offers.
-
BesselK.jl,
(MIT license) An ostensibly simple package that implements the modified
second-kind Bessel function $\mathcal{K}_\nu(x)$ natively in Julia for the
purpose of using automatic differentiation to compute derivatives with
respect to $\nu$. Which sounds very simple, but really isn't due to a
variety of subtle but legitimately problematic issues. See the paper
"Fitting Matern..." for details. The real point is to fit three-parameter
Matern covariance functions. The example file gives a simple Matern
covariance function that you can just drop in to your code (and your AD
stack, assuming you use ForwardDiff.jl) and then just forget about. But
please do cite the paper if you do that...