The principal aim of the paper is to investigate new integral expression $${\int}_0^{\infty}x^{s+1}e^{-{\sigma}x^2}L_m^{(\gamma,\delta)}\;({\zeta};{\sigma}x^2)\;L_n^{(\alpha,\beta)}\;({\xi};{\sigma}x^2)\;J_s\;(xy)\;dx$$, where $y$ is a positive real number; $\sigma$, $\zeta$ and $\xi$ are complex numbers with positive real parts; $s$, $\alpha$, $\beta$, $\gamma$ and $\delta$ are complex numbers whose real parts are greater than -1; $J_n(x)$ is Bessel function and $L_n^{(\alpha,\beta)}$ (${\gamma};x$) is generalized Laguerre polynomials. Some integral formulas have been obtained. The Maple implementation has also been examined.