Models where Stan outperforms Nutpie/Walnuts

I guess this is a bit of a detour from the main topic here, but let’s look at that function a bit more. I’m bad at naming things, but I’ll call it bisoftplus(x) = x + sqrt(1 + x ** 2) = exp(arcsinh(x)) for now. (because it behaves a bit like softplus for the inverse as well, so “bi”-directionally. Other ideas are welcome. Claude tells me someone proposed it as an activation function under the name of squareplus some time back). The inverse is (x^2 - 1) / 2x

We use maps from \mathbb{R} to \mathbb{R}^+ in two prominent places: As inverse link functions in regressions, and in the inverse unconstraining map. Those might look similar, but I think they often require different properties.

Constraining maps

The bisoftplus would be a pretty bad constraining map. If we wanted to sample a simple half-normal, the positive limit of the density as x goes to 0 turns into a cauchy-like left tail after an inverse bisoftplus transform. You can see how that happens if you think of the inverse of bisoftplus as sinh(log(x)). First, we take the log of the half-normal as we currently do, and then we exponentially blow up both tails with sinh. Not what we want.

softplus might have a case against exp in many applications. The left tail works like exp, but doesn’t squish the right tail of distributions. But the right tail of positive distributions is often perfectly fine as it is, there is no reason to squish it in exponentially.

But in constrast to softplus, the exp function does have some very friendly properties if we think of all of those functions as families of functions. Before we can log transform anything, we always have to remove the units, so there is a natural scaling factor we need to take into account. We also have the choice of which exponential base we want to use. That really gives us log(cx)/log(b) with the free hyperparameters c and b for the inverse-exp family. But those hyperparameters turn into a constant factor and a constant offset after the log transformation. HMC is invariant to an offset, and the diagonal mass matrix estimation will choose a scale for us anyway. So we can just use any values for c and b and completely ignore that those hyperparameters even exist. But this doesn’t work for softplus, where the scaling parameter really matters. We need to choose it somehow if we want the transformation to work properly. (We could actually learn it by minimizing the fisher divergence just as we learn the mass matrix/affine transform in nutpie, just food for thought…)

Link function for variance parameters

(cc @paul.buerkner from a previous discussion on this…)

So let’s say we want to run a regression, where we also predict variances based on some predictors. I’m guessing that this is roughly what’s happening in the model of @ssp3nc3r. The predictors live in \mathbb{R}, so we need an inverse link function to map that to the positive standard deviation.

I’ve seen exp and softplus in the wild for this, but they both have a lot of problems, and I think bisoftplus does something quite sensible here.

exp feels like a very natural choice, but it blows up predicted standard deviations. If our predictors have a normal distribution, then the standard deviation (and also the variance) has a log normal distribution. The log-normal has very long right tails, thicker than an exponential, but lighter than any pareto distribution (and often looks like a power law for quite a long time). This can lead to computational issues, but I don’t think it is usually what we want in a model either. “Sure, all standard deviations we’ve seen so far were between 0.5 and 2, but this out-of-sample standard deviation might just be 1000, who knows…”. I’ve never worked on anything where I thought long tails like this sound realistic.

The softplus fixes the long right tails, but it introduces new problems. We lose nice symmetry properties we had before. In the exp link, we had the property that if you swap the sign of a predictor, that the predicted standard deviation is 1/old. That’s gone. And the precision of the predicted values still has the very long tail of a log-normal, and that might be just as bad as the long tails of the standard deviation with the exp link.

The bisoftplus fixes both. We get our symmetry back, zero maps to 1, and the tails of the standard deviation and its inverse are friendly. And for large or small predictors it behaves additively on the standard deviation or its inverse. And around 0 it behaves like exp to second order.

One thing that wasn’t too important with the exp link is if we want to run the regression against the standard deviation or against the variance. That’s just a constant factor of 2 different with an exp link. But with a bisoftplus link it suddenly makes a difference. Do we want additivity of the standard deviation or the variance? Given that the variance is usually the thing that is additive, maybe it’s the better choice?

But I guess if we want to continue this, we should move it into a different thread.