First observe how the gamma distribution function can be written in terms of the incomplete gamma function.
$\boldsymbol{(1)} \qquad G(y) = \int_{0}^{y} \dfrac{c^{\gamma}}{\Gamma(\gamma)} x^{\gamma - 1}e^{-cx}dx = \dfrac{1}{\Gamma(\gamma)}\int_{0}^{cy} z^{\gamma- 1}e^{-z}dz = \mathfrak{S}(\gamma,cy )$
$Y \sim \Gamma(\gamma,c)$
Now, in my book the statement is: "From this we see that the family of gamma distributions is closed towards multiplication with a positive constant, that is, for $\rho > 0 $ we have $\rho Y \sim \Gamma(\gamma,c/\rho)$"
My question is: How can i see from (1) that the gamma distribution has this property (closed toward multiplication)
You can do it directly: let $W = \rho Y$ where $Y$ is distributed as described. Then $$\Pr[W \le w] = \Pr[\rho Y \le w] = \Pr[Y \le w/\rho] = G(w/\rho).$$ Then we have $$G(w/\rho) = \mathfrak{S}(\gamma,cw/\rho) = \mathfrak{S}(\gamma,(c/\rho)w).$$ Thus $W \sim \Gamma(\gamma,c/\rho)$. All of this is just a fancy way of showing that $c$ is a scale parameter of the gamma distribution; that is to say, if $Y \sim \operatorname{Gamma}(\gamma,c)$, then $\rho Y \sim \operatorname{Gamma}(\gamma,c/\rho)$, if using the rate parametrization for the gamma distribution.