How many real solutions does $\lambda_1 e^{y - \lambda_1 e^y} (1 - \lambda_1 e^y ) + \lambda_2 e^{y - \lambda_2 e^y} (1 - \lambda_2 e^y ) = 0$ have? I am trying to work out for what $\lambda_1, \lambda_2 > 0$ is it true that  $f(y) = \lambda_1 e^{y-\lambda_1 e^y} + \lambda_2 e^{y-\lambda_2 e^y}$ is unimodal?
Experimentally it seems it is unimodal when $\lambda_1 < \lambda_2$ and $\frac{\lambda_2}{\lambda{1}} <  7.5$ .
To work this out I started with:
$$\frac{d}{dy} \left(\lambda_1 e^{y-\lambda_1 e^y} + \lambda_2 e^{y-\lambda_2 e^y} \right) = \lambda_1 e^{y - \lambda_1 e^y}  (1 - \lambda_1 e^y ) + \lambda_2 e^{y - \lambda_2 e^y}  (1 - \lambda_2 e^y )$$
It seems we then need to check when
$$\lambda_1 e^{y - \lambda_1 e^y}  (1 - \lambda_1 e^y ) + \lambda_2 e^{y - \lambda_2 e^y}  (1 - \lambda_2 e^y ) = 0$$
has more than one solution when solved for $y \in \mathbb{R}$.  How can we determine the conditions under which it has different numbers of solutions?
Added:
Substituting $z = e^y$ and dividing by $e^{y-1}$ we are trying to determine how many solutions
$$
 \lambda_1 e^{1-\lambda_1 z}(1-\lambda_1 z) +\lambda_2e^{1-\lambda_2 z}(1-\lambda_2 z) = 0
$$
has with $z > 0$.
Examples:
Example $\lambda_1 = 1, \lambda_2 = 7$ with only one mode (code in python):
import matplotlib.pyplot as plt
import numpy as np
def pdf_func(y, params):
    return sum([lambd*np.exp(y - lambd * np.exp(y)) for lambd in params])
params = [1, 7]
xs = np.linspace(-10,10,1000)
plt.plot(xs, [pdf_func(y, params) for y in xs])


Example $\lambda_1 = 1, \lambda_2 = 50$ with two modes:

Questions

*

*How can one prove (assuming it is true) that that the number of local maxima that $f(y)$ has is either 1 or 2 and there are no other possibilities?

*Is it true that for $\lambda_2 > \lambda_1 > 0$, there exists a threshold $c$ so that if $\frac{\lambda_2}{\lambda_1} < c$ then $f(y)$ is unimodal and if not it has two local maxima? (My guess is that the answer is yes and this threshold is around $7.5$.)

 A: Define $f(x)=ae^{-ax}(1-ax)+be^{-bx}(1-bx)$ for $x>0$, where $b>a>0$. 
The OP is interested in the roots of $f$.
Let $r=\frac ba$.

Lemma 1 $f$ has exactly one root in $(0,\frac2b)$.

Sketch of proof:


*

*$f(0)=a+b >0$

*$\displaystyle{f\left(\frac2b\right)<0}$ 

*The derivative of $g_c(x):= ce^{-cx}(1-cx)$ is $g_c’(x)=-c^2e^{-cx}(2-cx)$. Thus, $g_c$ is decreasing on $(0,\frac2c)$, and $f=g_a+g_b$ is decreasing on $(0,\min(\frac2a,\frac2b))=(0,\frac2b)$.


By intermediate value theorem, there exists at least one zero in $(0,\frac2b)$. Monotonicity implies uniqueness of the zero.

Lemma 2
$$\text{number of zeroes in $\left[\frac2b,\infty\right)$}=
\begin{cases}
0, & r<\rho \\
1, & r=\rho \\
2, & r>\rho \\
\end{cases}
$$
  where $\rho\approx 7.566$ and satisfies $$\exp\left[h(\rho)\left(1-\frac1\rho\right)\right]+\rho^2\cdot\frac{1-h(\rho)}{\rho-h(\rho)}=0,\quad h(t)=\frac12\left(1+t+\sqrt{t^2-6t+1}\right)$$

Ideas:
I have not been able to come up with a formal proof. But intuition and numerical experiments have led me to believe that the critical case occurs when the zero of $f$ is also a stationary point.



Therefore, we want to solve for $\rho=\frac ab$, given that
$$\quad f(x)=ae^{-ax}(1-ax)+be^{-bx}(1-bx)=0$$
$$\implies ae^{-ax}(1-ax)=-be^{-bx}(1-bx)\qquad{(*)}$$
$$\quad f'(x)=-a^2e^{-ax}(2-ax)-b^2e^{-bx}(2-bx)=0$$
$$\implies a^2e^{-ax}(2-ax)=-b^2e^{-bx}(2-bx)$$
Dividing the two equations give
$$a\cdot\frac{2-ax}{1-ax}=b\cdot\frac{2-bx}{1-bx}$$
With a little algebra we get
$$ab\cdot x^2-(a+b)x+2=0$$
By quadratic formula,
$$x=\frac12\left(\frac1a+\frac1b+\sqrt{\frac1{a^2}+\frac1{b^2}-\frac6{ab}}\right)
=\frac1{2b}\left(1+\rho+\sqrt{\rho^2-6\rho+1}\right)=\frac{h(\rho)}b$$
(Again, experiments tell that only the positive root should be considered.)
Hence, $ax=h(\rho)/\rho$ and $bx=h(\rho)$. Substituting this back into $(*)$,
$$ae^{-h(\rho)/\rho}\left(1-\frac{h(\rho)}{\rho}\right)=-be^{-h(\rho)}(1-h(\rho))$$
or
$$\exp\left[h(\rho)\left(1-\frac1\rho\right)\right]+\rho^2\cdot\frac{1-h(\rho)}{\rho-h(\rho)}=0$$
This equation has two roots, $\rho\approx7.566$ or $\rho’=\frac1\rho\approx 0.132$. The latter is rejected as $b>a\implies \rho>1$.
Any ideas on how to turn these observations into a formal proof?

23rd October 2019 edit
I'd like to view the whole problem from a slightly different perspective. Alternatively, consider the system
$$
\begin{cases}
ye^{-yx}(1-yx)+ae^{-ax}(1-ax)=0 \\
y=b \\
y\ge a
\end{cases}
\qquad{(S)}
$$
One advantage of this approach is that the increasing trend of the number of solution ($1\to 2\to 3$) is more obviously illustrated by the implicit function:

The proof may be thus simplified.
Intuition: What does this graph actually mean? Imagine that you are plotting $y=be^{-bx}(1-bx)+ae^{-ax}(1-ax)$ for a fixed $a$ and a varying $b$. Then:


*

*Take a snapshot of the x-axis of the graph.

*Make a small increment to $b$.

*Take another snapshot of the x-axis of the graph.

*Stack/append this snapshot on top of the previous one.

*Repeat.


Then you will get the graph above.
Since we are interested only in the distribution of x-intercepts (i.e. roots) of $y=be^{-bx}(1-bx)+ae^{-ax}(1-ax)$ for different values of $b$, the graph above essentially captures all the information that is of our interest.
A few more words: The implicit function has two branches, and they are separated by $x=\frac2y$. (This observation is indeed a direct consequence of Lemma 1.)
In other words, 
$$
\begin{cases}
ye^{-yx}(1-yx)+ae^{-ax}(1-ax)=0 \\
y>\frac2x \\
y\ge a
\end{cases}
$$
uniquely defines a 'function', which is the parabola-like branch in the graph. 
We will focus on this branch only. Call this branch $Y_a(x)$.

To prove Lemma 2, it suffices to prove the equivalent version

$$\text{The number of solutions of $Y_a(x)=b$ is }
\begin{cases}
0, & r<\rho \\
1, & r=\rho \\
2, & r>\rho \\
\end{cases}
$$

The basic steps of the proof will be:


*

*Prove $\lim_{x\to 0^+}Y_a(x)=+\infty$.

*Prove $\lim_{x\to (1/a)^-}Y_a(x)=+\infty$.

*$Y_a(x)$ has only one stationary point in $(0,\frac1a)$, and it is a minimum attained at $y=a\rho$.


When these three statements are proved, it can be shown that


*

*$Y_a(x)=b$ has no solutions when $r<\rho$, as $b$ is smaller than the minimum of $Y_a$.

*$Y_a(x)=b$ has exactly one solution when $r=\rho$, and the solution is at the minimum.

*$Y_a(x)=b$ has two solutions when $r>\rho$. Let the coordinates of the minimum be $(k,a\rho)$. Consider $Z(x)=Y_a(x)-b$


*

*Since $Z(0)=+\infty>0$ and $Z(k)=a\rho-b<0$, by intermediate value theorem (IVM) $Z(x)$ has roots in $(0,k)$. Moreover, absence of local minimum in $(0,k)$ implies monotonic decrease of $Z$. Therefore, $Z(x)$ has one unqiue root in $(0,k)$.

*Similarly, $Z(k)<0$ and $Z(\frac1a)=+\infty>0$ together with absence of local maximum in $(k,\frac1a)$ give another unqiue root in $(k,\frac1a)$.



This is the full outline of the proof. What remains is the proof of the first three statements.

Here is the proof of the first statement:
From the definition of $Y$, $$Y>\frac2x\implies\lim_{x\to 0^+}Y>\lim_{x\to 0^+}\frac2x=+\infty$$
Therefore, $Y$ diverges to $+\infty$ as $x\to0^+$.
Proof of statement 2:
Let $\phi(z)=ze^{-z}(1-z)$.
Then $\phi(xY)+\phi(ax)=0$.
Assume $Y$ is bounded in a neighbourhood of $x=\frac1a^-$.
Taking limits on both sides,
$$\lim_{x\to(1/a)^-}\phi(xY)=0$$
$$\phi\left(\frac1a\lim_{x\to(1/a)^-}Y\right)=0$$
$$\implies Y\to0\text{ or }Y\to a$$
Either value does not satisfy the definition of $Y$: $Y>\frac2x$.
Hence $Y$ is not bounded near $x=\frac1a^-$. If $Y$ diverges to $-\infty$, then $\lim_{x\to(1/a)^-}\phi(xY)=-\infty\ne 0$. Hence, $Y$ diverges to $+\infty$.
The proof of statement 3 will be added soon.
A: Consider 
$$
f(z) = \lambda _1 z e^{-\lambda _1 z}+\lambda _2 z e^{-\lambda _2 z}
$$
with $z = e^y$
and now
$$
f'(z) = -\lambda _1 e^{-\lambda _1 z} \left(\lambda _1 z-1\right)-\lambda _2 e^{-\lambda _2 z} \left(\lambda
   _2 z-1\right)=0
$$
or
$$
\frac{\lambda_2^2}{\lambda_1^2}e^{(\lambda_1-\lambda_2)z}+\frac{z-\frac{1}{\lambda_2}}{z-\frac{1}{\lambda_1}} = 0
$$
here
$$
\frac{\lambda_2^2}{\lambda_1^2}e^{(\lambda_1-\lambda_2)z}\gt 0
$$
so the zero's location is associated to the sign of
$$
\frac{z-\frac{1}{\lambda_2}}{z-\frac{1}{\lambda_1}}
$$
so
$$
\min_i\frac{1}{\lambda_i}\le z \le \max_i\frac{1}{\lambda_i}
$$
NOTE
As an example, for $\lambda_1 = 1, \lambda_2 = 50$ we have the sign change plot

Now the solutions should be searched into the negative interval which is $0.02\le z \le 1$. The roots here are $z =\{0.0211011,0.113856,1\}$
Now considering after $\lambda_2 = \zeta\lambda_1$
$f'(z,\zeta,\lambda_1)=\lambda _1 \zeta  e^{-\lambda _1 z \zeta } \left(\lambda _1 z \zeta -1\right)+\lambda _1 e^{-\lambda _1 z}
   \left(\lambda _1 z-1\right)$. 
or calling $y = \lambda_1 z, x = \zeta$
$$
g(x,y) = x(x y -1)e^{-y(x-1)}+y - 1
$$
or also
$$
x y +W\left(\frac{(y-1)e^{1-y}}{x}\right) - 1 = 0
$$
Here $W(\cdot)$ is the Lambert function.
Follows the plot for $g(x,y) = 0$

As can be observed roughly for $0\le x \le 7.5$ we have one root and for $7.5 \lt x$ we have three roots.
The determination of the passage from one root to three is made as follows. Solving numerically $g(x,y) = 0$ and $\frac{dg}{dy} = 0$ we obtain the point $(7.56628, 0.802977)$ (the black point as the intersection of $g(x,y) = 0$ in blue and $\frac{dg}{dy} = 0$ in green)

