When plotting over a big range, especially when plotting power laws or things with exponential dependence, the logarithm function is your friend. The log-log plot and semi-log plot are standard tools for visualization in science.

But while the logarithm function works on the domain (0,∞), many applications in science operate on other domains, like (-∞,∞) or (0,1).  In these cases, logarithmic plotting is not always appropriate, but in some cases it can be and it can be frustrating to find a good visualization tool that captures what you want.

For instance, I’m teaching stellar structure and evolution again and it’s tricky to make a plot that captures the details of the stellar atmosphere, which is only 1 part in 10 billion of the mass, and the details of the core,  where you want to see details covering only, say, the last 0.1% of the mass.

Let’s say you wanted to compare, for instance, the radiative temperature gradient ∇_rad in the sun to the actual temperature gradient (this difference is related to the efficiency of convection). This quantity spans a huge range in the very thin outer layer of the star, so trying to plot its dependence as a function of mass inside the sun gives a pretty useless result:

Maybe we just need to plot logarithmically?  Let’s try it:

That definitely helps!  We see that there is stuff going on down in the core, and apparently the gradient returns to around 1 near the surface, but we still can’t really tell what’s going on.

There are lots of solutions—you can plot in radius instead of mass, or use log(Pressure), for instance, or perhaps we could take the log of 1-m/M?  That would expand the region of interest logarithmically.  Let’s try that:

Now we’re talking! Suddenly we can see all of the action in the convective envelope.  But…we’ve lost resolution on the core.  Is there some way to expand both ends logarithmically?

There certainly is! We can use the logit function, which I’ve written about it before as a way to expand things on the domain (0,1). This is just the sum of log(m/M) and -log(1-m/M). [nb: logit is usually defined base-e, so you’ll have to divide it by ln(10) to get logit₁₀].

This has the lovely property that is logarithmic in both directions, near zero and near 1.  So a value of -5 means I’m at 10^-5, and a value of 5 means I’m at 1-10^-5, or 0.99999.  In other words, negative vales “count zeroes”, and positive values “count nines”.  0 corresponds to a half, right in the middle.  So it’s intuitive to read and understand and gives you the dynamic range you want at both ends.  It’s a great way to plot things like transmission or conductivity when you have values that hug 0 and 1, and you care about how close they get.

Let’s try it:

Perfect!  Now the core is expanded way way out so we can see what’s going on (for the experts: it’s slightly convective here because this is the ZAMS sun and we have extra CNO fusion going on until all the carbon burns out).

There are other situations where you have functions that have power law dependence in both the negative and positive directions, as well. This has come up occasionally with me, and I’m always tempted to do something like sign(x)log₁₀(|x|), but it doesn’t work because near zero it blows up.

What you want is a function that is linear near zero, and asymptotes to sign(x)log₁₀(|x|).

It turns out matplotlib has a ‘symlog’ scaling that addresses this!  I’m not totally sure how it works, but here it is. This is a linear-symlog plot of a function against itself (so, the 1:1 line):

This is nice, but there’s a clear kink in the data where it switches from logarithmic to linear in between.  It’s pretty kludgy.

But it turns out that there is an already-defined not-at-all kludgy function that already does this smoothly!  It’s the inverse hyperbolic sine function!  Specifically, in base-10 you’d call arcsinh(x/2)/ln(10), which is also log₁₀(x/2+√(x/2)²+1)). And it behaves just like you want:

BTW, it’s not just that arcsinh has a sigmoid-y kind of shape so gets the job done, it’s actually the optimal function here.  That’s because it is the inverse of sinh, which is exactly the function that should generate a straight line in this scaling: sinh(x) is the average of exp(x) and -exp(-x), which are the functions you want on either side of zero; far from zero the other one decays away, and at zero the function is nicely linear (the second order term is zero).

[Update: Prior art from Warrick Ball, who implies this is on old trick and super useful for color scales!:

https://twitter.com/warrickball/status/1040199586344783879?s=20

]

I don’t know what to call this scaling (inverse hyperbolic sine is too long and technical and obscure for what it does) but I’ll suggest bilogarithmic which seems to be a term that people sometimes use to mean “logarithmic in both directions” but has no well established meaning.

Anyway that’s how to get functions that are logarithmic in both directions, for domains that are bounded on both, one, or no sides!

The log scaling is already in matplotlib, of coure. There is also a ‘logit’ scaling but it does not work very well for me in the above examples, so could use some tweaking—in particular it would be good to have a base-10 version.  And ‘symlog’ could be replaced with arcsinh(x)/log(10).  This would reduce some of its functionality (right now you choose the range over which to make it linear) but I think it’s worth it to have it better behaved at the kink.

I don’t know enough Python to implement them myself, but perhaps it can be a feature request if enough people start using these!