I know that the value of the integral is as follows $$\int_0^{\infty}\frac{(1-e^{-\lambda z})}{\lambda^{a+1}} d \lambda =z^a \frac{\Gamma(1-a)}{a}$$
However, how exactly the integral is calculated? How one proves that the integral even exists? Why does the integral only exist for $0<a<1$?
Existence: at $\lambda=\infty$, $a>0$ guarantees that the integrand is below $1/\lambda{1+a}$ which is integrable. At $\lambda=0$, $$ \int_0^\delta\frac{(1-e^{-\lambda z})}{\lambda^{a+1}}\,d\lambda=\int_0^\delta\frac{\lambda z+O(\lambda^2)}{\lambda^{a+1}}\,d\lambda=\int_0^\delta\frac{z}{\lambda^a}+O(\lambda^{1-a})\,d\lambda. $$ The first term requirest $a<1$ for convergence, and the second one $a<2$. So $0<a<1$ are precisely the values of $a$ for which your integral exists.
Regarding the calculation, I wouldn't know how to find the formula, but here is a way to justify it. It suffices to use differentiation under the integral sign: we have $$ \frac d{dz}\,\int_0^\infty\frac{(1-e^{-\lambda z})}{\lambda^{a+1}}\,d\lambda=\int_0^\infty\frac{\lambda\,e^{-\lambda z}}{\lambda^{a+1}}\,d\lambda=\int_0^\infty\lambda^{-a}\,e^{-\lambda z}\,d\lambda=\int_0^\infty z^{a-1}\,t^{-a}\,e^{-t}\,dt=z^{a-1}\Gamma(1-a)=\frac d{dz}\,\left(\frac{z^a}a\,\Gamma(1-a)\right). $$ So the two sides of your equality differ by a constant $c$. When $z=0$ both sides equal to $0$, so $c=0$.