Sum of exponential random variable with different means

$\newcommand{\+}{^{\dagger}}% \newcommand{\angles}[1]{\left\langle #1 \right\rangle}% \newcommand{\braces}[1]{\left\lbrace #1 \right\rbrace}% \newcommand{\bracks}[1]{\left\lbrack #1 \right\rbrack}% \newcommand{\ceil}[1]{\,\left\lceil #1 \right\rceil\,}% \newcommand{\dd}{{\rm d}}% \newcommand{\ds}[1]{\displaystyle{#1}}% \newcommand{\equalby}[1]{{#1 \atop {= \atop \vphantom{\huge A}}}}% \newcommand{\expo}[1]{\,{\rm e}^{#1}\,}% \newcommand{\fermi}{\,{\rm f}}% \newcommand{\floor}[1]{\,\left\lfloor #1 \right\rfloor\,}% \newcommand{\half}{{1 \over 2}}% \newcommand{\ic}{{\rm i}}% \newcommand{\iff}{\Longleftrightarrow} \newcommand{\imp}{\Longrightarrow}% \newcommand{\isdiv}{\,\left.\right\vert\,}% \newcommand{\ket}[1]{\left\vert #1\right\rangle}% \newcommand{\ol}[1]{\overline{#1}}% \newcommand{\pars}[1]{\left( #1 \right)}% \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\pp}{{\cal P}}% \newcommand{\root}[2][]{\,\sqrt[#1]{\,#2\,}\,}% \newcommand{\sech}{\,{\rm sech}}% \newcommand{\sgn}{\,{\rm sgn}}% \newcommand{\totald}[3][]{\frac{{\rm d}^{#1} #2}{{\rm d} #3^{#1}}} \newcommand{\ul}[1]{\underline{#1}}% \newcommand{\verts}[1]{\left\vert\, #1 \,\right\vert}$ With $\large\ t > 0$: \begin{align} &\color{0000ff}{\large{\rm P}\braces{x + y < t}} \\[3mm]&=\int_{0}^{\infty}\lambda\expo{-\lambda x}\int_{0}^{\infty} \mu\expo{-\mu y}\Theta\pars{t - x - y}\,\dd x\,\dd y = \lambda\mu\int_{0}^{\infty}\expo{-\lambda x}\int_{0}^{\infty} \expo{-\mu y}\Theta\pars{\bracks{t - x} - y}\,\dd y\,\dd x \\[3mm]&= \lambda\mu\int_{0}^{\infty}\expo{-\lambda x}\Theta\pars{t - x}\int^{t - x}_{0} \expo{-\mu y}\,\dd y\,\dd x = \lambda\mu\int_{0}^{\infty}\expo{-\lambda x}\Theta\pars{t - x}\, {\expo{-\mu\pars{t - x}} - 1 \over -\mu}\,\dd x \\[3mm]&= \lambda\int_{0}^{t} \bracks{\expo{-\lambda x} - \expo{-\mu t}\expo{\pars{\mu - \lambda}x}}\,\dd x = \lambda\bracks{% {\expo{-\lambda t} - 1 \over -\lambda} - \expo{-\mu t}{\expo{\pars{\mu -\lambda}t} - 1 \over \mu - \lambda}} \\[3mm]&= \lambda\bracks{% {1 - \expo{-\lambda t} \over \lambda} + {\expo{-\mu t} - \expo{-\lambda t} \over \mu - \lambda}} = \lambda\, {\pars{\mu - \lambda}\pars{1 - \expo{-\lambda t}} + \lambda\pars{\expo{-\mu t} - \expo{-\lambda t}} \over \lambda\pars{\mu -\lambda }} \\[3mm]&= {\pars{\mu - \lambda} - \pars{\mu - \lambda}\expo{-\lambda t} + \lambda\pars{\expo{-\mu t} - \expo{-\lambda t}} \over \mu -\lambda} =\color{#0000ff}{\large% 1 - {\mu\expo{-\lambda t} -\lambda\expo{-\mu t} \over \mu - \lambda}} \end{align}

Felix has provided a very nice detailed solution by hand. This type of problem can also be easily solved by automated methods using a computer algebra system ...

By independence, the joint pdf of $(X,Y)$ is say $f(x,y)$:

We seek the cdf of the sum, i.e. $P(X+Y<t)$:

where I am using the Prob function from the mathStatica add-on to Mathematica to calculate the probability automatically. All done.

Notes:

1. The answer is the same as that obtained by Felix.
2. As formal disclosure, I should perhaps add that I am one of the authors.