How do I smoothly merge two power laws? I've got data that, when plotted on a log-log scale, looks like two straight lines that gently bend and blend into each other. See the picture below.
I don't know how to find a function to fit that. Now, I know that I should use some kind of broken power law, and I know how to make one that asymptotically converges to a constant, but I don't know how to combine two power laws. I've found something that looks just right as far as its shape is concerned, and even some Python code to generate it, but I can't figure out how to tweak it to fit my data without breaking the continuity between the two power laws:
def chabrier03individual(m):
        k = 0.158 * exp(-(-log(0.08))**2/(2 * 0.69**2))
        return numpy.where(m <= 1,
                0.158*(1./m) * exp(-(log(m)-log(0.08))**2/(2 * 0.69**2)), k*m**-2.3)

Note that here, log means $log_{10}$. Am I even on the right track?

 A: I've always liked the very simple
$$y=\frac{1}{\frac{1}{x^a}+\frac{1}{x^{-b}}}$$
where both $a$ and $b$ are positive.
When $x$ is small, the $x^a$-term dominates and you get a slope of $a$ on the log-log plot, whereas when $x$ is large the $x^{-b}$-term dominates with a slope of $-b$ on the log-log plot.
Of course this can be modified with additional constant multipliers in the denominator.
A: After many years this question was posted, I came across this problem and found this thread. However, my problem is more general: how do I smoothly merge two or more power laws?
The goal is to find a function that has a derivative ($\frac{d\log y}{d\log x}$) that looks like a steps function with the discontinuity smoothly joined. The most natural candidate is the logistic function:
$$
f(x) = \frac{L}{1 + \exp(-k(x - x_0))}.
$$
Given sequences of power-law indices $m_1, m_2, ..., m_n$ and turn-over points (in ascending order) $x_{1,2}, x_{2,3}, ..., x_{n-1,n}$, choosing $k = 1$, we can write the derivative as
$$
\frac{d\log y}{d\log x} = m_1 + \frac{m_2 - m_1}{1 + \exp(-\log x + \log x_{1,2})} + ... + \frac{m_n - m_{n-1}}{1 + \exp(-\log x + \log x_{n-1,n})}.
$$
It turns out that this differential equation can be rewritten to a very simple form:
$$
\frac{d\log y}{dx} = \frac{m_1}{x} + \frac{m_2 - m_1}{x + x_{1,2}} + ... + \frac{m_n - m_{n-1}}{x + x_{n-1,n}}.
$$
After integration, it reduces to a very simple form that one should have thought of at the beginning:
$$
y = A x^{m_1} (x + x_{1,2})^{m_2 - m_1} ... (x + x_{n-1,n})^{m_n - m_{n-1}}
$$
The reason it works is that, for some $x$ that falls between two turn-over points $x_{i-1,i}$ and $x_{i,i+1}$, $x$ is much larger than $x_{j-1,j}$ for all $j \leq i$ while $x$ is much smaller than $x_{j,j+1}$ for all $j \geq i$. Therefore, all former terms reduces to $x^{m_i}$ and all latter terms reduces to a constant. And there you have it!
A: What you provided in the previous question looked like a the cumulative data function (CDF), or survival function (SF), and this the probability density function (PDF). So this is a broken power law with a single break. You can either use the negative derivative of the SF from the previous question, which gives:
\begin{align}
f(x) & \equiv - N\frac{\partial B(x, s, \beta, x_n)}{\partial x} \\
& = - N\beta \left[1 + \left(\frac{x}{x_n}\right)^{|\beta|s}\right]^{\operatorname{sgn}(\beta)/s - 1} \left(\frac{x^{|\beta|s-1}}{x_n^{|\beta|s}}\right)^{|\beta|s},
\end{align}
or, if you'd prefer to fit the pdf with its own smoothly broken power law that has a single break, you'd go for:
$$
f(x) = N x^n \left[1+\left(\frac{x}{x_n}\right)^{|\beta|s}\right]^{\operatorname{sgn}(\beta)/s},
$$
which gives one more parameter than fitting the SF with a broken power law.
