Some explanations of the sequence of derivations\begin{align}
P(Y\le y) &\stackrel{\text{(1)}}{=} P(Y_1 + Y_2<y)\\
&\stackrel{\text{(2)}}{=} \mathbb E[\mathbb I_{Y_1 + Y_2<y}]\\
&\stackrel{\text{(3)}}{=} \mathbb E[\mathbb E[\mathbb I_{Y_1 + Y_2<y}|Y_2]]\\
&\stackrel{\text{(4)}}{=} \int_0^{\infty} \mathbb E[\mathbb I_{Y_1 + Y_2<y}|Y_2=u]\,f(u;\theta)\,\text du \\
&\stackrel{\text{(5)}}{=} \int_0^{y} \mathbb E[\mathbb I_{Y_1 + Y_2<y}|Y_2=u]\,f(u;\theta)\,\text du \\
&\stackrel{\text{(2)}}{=} \int_0^{y} P(Y_1 + Y_2 \le y |Y_2 = u)f(u;\theta)\,\text du \\
&= \int_0^{y} P(Y_1 + Y_2 \le y |Y_2 = u)\frac{1}{\theta}e^{\frac{-u}{\theta}}\,\text du\\
&\stackrel{\text{(6)}}{=} \int_0^y P(Y_1 \le y -u |Y_2 = u)\frac{1}{\theta}e^{\frac{-u}{\theta}}\,\text du \\
&\stackrel{\text{(7)}}{=} \int_0^y P(Y_1 \le y -u)\frac{1}{\theta}e^{\frac{-u}{\theta}}\,\text du \\
&\stackrel{\text{(8)}}{=} \int_0^y \left(1 - e^{\frac{-(y-u)}{\theta}}\right)\frac{1}{\theta}e^{\frac{-u}{\theta}}\,\text du \\
&= \int_0^y \frac{1}{\theta}e^{\frac{-u}{\theta}}\,\text du - \int_0^y e^{\frac{-(y-u)}{\theta}}\,\frac{1}{\theta}e^{\frac{-u}{\theta}}\,\text du \\
&= P(Y_2\le y) - \int_0^y \frac{1}{\theta}e^{\frac{-(y-u+u)}{\theta}}\,\text du\\
&\stackrel{\text{(9)}}{=} 1 - e^{\frac{-y}{\theta}}-\frac{1}{\theta}ye^{\frac{-y}{\theta}}\end{align}
where
- by definition of $Y=Y_1+Y_2$
- when writing the probability of an event as the expectation of the indicator of this event
- when using the law of total expectation
- by definition of the expectation
- when accounting for the fact that both $Y_1$ and $Y_2$ are positive rv's
- because $Y_2$ is known to be equal to $y$
- because $Y_1$ and $Y_2$ are independent $\mathcal E xp(1/\theta)$
- as there was a mistake of an extra $1/\theta$ in the corresponding row from the question (before my editing), since the cdf of $Y_1$ is $F(y)=1-\exp\{-y/\theta\}$
- the final result is nonetheless correct.
An alternative approach is to determine first the density of $Y=Y_1+Y_2$ using the convolution theorem: the density of $Y$ writes as
\begin{align}
f_Y(y) &= \int f_{Y_1}(u)f_{Y_2}(y-u)\,\text du\\
&= \int_0^\infty f_{Y_1}(u)f_{Y_1}(y-u)\,\text du\\
&= \int_0^y \frac{1}{\theta}e^{-u/\theta}\,\frac{1}{\theta}e^{-(y-u)/\theta}\,\text du\,\mathbb I_{y>0}\\
&= \frac{1}{\theta^2} \int_0^y e^{-y/\theta}\,\text du\,\mathbb I_{y>0}\\
&= \frac{ye^{-y/\theta}}{\theta^2}\,\mathbb I_{y>0}
\end{align}
which is a Gamma$(2,1/\theta)$ density. The probability (cdf) then follows
\begin{align}
P(Y<y)&=\int_0^y \frac{ue^{-u/\theta}}{\theta^2}\,\text du\\
&=\int_0^y -\frac{u}{\theta}\frac{\text de^{-u/\theta}}{\text du}\,\text du\\
&=-\frac{y}{\theta}e^{-y/\theta}+\int_0^y \frac{\text du}{\text du}
\frac{e^{-u/\theta}}{\theta}\,\text du\\
&= -\frac{1}{\theta}ye^{\frac{-y}{\theta}} + 1 - e^{\frac{-y}{\theta}}
\end{align}