The data is showing power law behaviour, and although empirical power laws are often poorly understood, in this case there is a plausible rough explanation based on theoretical relationships between decay energy and lifetime for alpha and beta decays.
Here are NuDat 3.0 data for 1797 radionuclides plotted in the same manner as in the question. (I get a 404 error when I try to access the data using the question's repository code.)
This plot is very similar to the plot in the question, but the two steps at long lifetimes in that plot do not appear and are almost certainly artifacts. It is hard to say more about the steps without the question's actual source data.
Instead of a cumulative plot, this probability density plot of the same data better shows the behaviour we are trying to understand:
The fuzzy jaggedness of the data (shown in grey) is because the empirical probability density is calculated directly from the data and is not smoothed or fit.
The distribution falls off roughly as as the inverse of the half-life, i.e. $t_{h}^{-1}$.
A theoretical explanation of why we would roughly expect a $1/t_h$ can be found in Corral, Font, and Camacho's paper on "Non-characteristic Half-lives in Radioactive Decay".
Alpha and beta decays need to be considered separately.
Alpha decay is a quantum tunnelling process, and the decay energy $Q$ and half-life $t_{h}$ are related by the Geiger-Nuttal Rule, which we can write as:
$$
\ln{\frac{t_{h}}{A}} = \frac{BZ}{\sqrt{Q}} = \frac{B}{\sqrt{U}}
$$
where $Z$ is the radionuclide's atomic number, $Q$ is the decay's energy, $U\equiv Q/Z^2=B^2/(\ln^2 t_h/A)$, and $A$ and $B$ are coefficients that are approximately constant. The (unnormalized) probability density can then be written
$$
D_\alpha = \frac{dN_{\alpha}}{dU}\left|{\frac{dU}{dt}}\right|=f\left(\frac{B^2}{\ln^2 t_h/A}\right)\left(\frac{2B^2}{\ln^3 t_h/A}\right)\frac{1}{t_h}
$$
We don't know the form of $f$, but because it only depends logarithmically on half-life, it is a slowly varying function for pretty much any reasonable choice. The only non-logarithmic half-life dependence in the expression is the final $1/t_h$ factor, so this is expected to dominate.
Beta decay is a four-fermion interaction that roughly follows Sargent's Rule:
$$
t_{1/2} \sim \frac{C^5}{Q^5}
$$
and the probability density can be written as
$$
D_\beta \sim \frac{dN_{\beta}}{dQ}\left|{\frac{dQ}{dt}}\right|=g\left(\frac{C}{t^{1/5}}\right)\frac{1}{t_h^{1+1/5}}
$$
If $g$ is a slowly varying function of $t_h$, then this will give a $1/t^{1.2}$ power law, but the validity of the "slowly varying" assumption is less obvious because the coefficient $C$ varies by orders of magnitude for different radionuclides. According to the Corral-Font-Comacho paper, however, the Pickands–Balkema–De Haan statistical theorem implies that $g$ is likely to have a Pareto power-law tail, and the combined result is that $D_\beta$ is likely to have roughly inverse power law behaviour with an exponent between $1$ and $1.2$.
So it is reasonable that the observed distribution has a power law distribution with exponent around $1$. As shown in the figure, the exponent seems to vary from $\sim 0.9$ to $\sim 1.2$. It is not immediately obvious if the exponent $<1$ for short half-lives is due to physics or some sort of selection bias.
