How to rewrite this integral as some non-elementary function or any possible close-form $?$: $$ I \equiv \int_{0}^{p}2\left(a \over y\right)^{M/2}\ \operatorname{K}_{M}\left(2\,\sqrt{\,{b \over y}}\,\right) \exp\left(-\,{y \over c}\right)\,\mathrm{d}y $$ where $a,b,c,p$ are all positive integer and $\operatorname{K}_{M}\left( {} \right)$ is the modified Bessel function of the second kind with $M$ is a positive integer that can be $1,2\ \mbox{or}\ 3$.
Currently, I am not considering the general case for all positive integer $M$ but any contribution is great.
The closest identity I could find that might be useful is:
$$\int_0^\infty t^{\mu-1}~e^{-at}~K_{\nu}(t)~\mathrm{d}t=\begin{cases} \sqrt{\frac{\pi }{2}} \Gamma ( \mu -\nu ) \Gamma ( \mu +\nu )\left( 1-a^{2}\right)^{-\mu /2+1/4} P^{-\mu +1/2}_{-\nu -1/2}( a) & -1< a< 1\\ \sqrt{\frac{\pi }{2}} \Gamma ( \mu -\nu ) \Gamma ( \mu +\nu )\left( 1-a^{2}\right)^{-\mu /2+1/4} P^{-\mu +1/2}_{\nu -1/2}( a) & \operatorname{Re}( a) \geq 0\ ;\ a\neq 1 \end{cases}$$
With $P^\alpha_\beta$ the associated Legendre function.
REFERENCE: https://dlmf.nist.gov/10.43