For arbitrary real a>0 this is a special case of the generalized _p\Psi_q(A;B;ζ) Fox-Wright function, where A=[(a_1,\alpha_1),(a_2,\alpha_2),...,(a_p,\alpha_p)] and B=[(b_1,\beta_1),(b_2,\beta_2),...,(b_q,\beta_q)] being a_j, j=1,..,p and b_k, k=1,..,q complex parameters and \alpha_j, \beta_k are positive. _p\Psi_q(A;B;ζ)=\sum_{n=0}^\infty \frac{ζ^n}{n!}\frac{\prod_{j=1}^{p}\Gamma(a_j+\alpha_jn)}{\prod_{k=1}^{q}\Gamma(b_k+\beta_kn)} None gamma function in the numerator is singular. This means a_{j}+α_{j}m≠-ℓ, with j=1,2,..,p ∧ ℓ,m∈ℕ₀. Series convergence depends on \kappa, \rho, ϑ κ=∑_{j=1}^{q}β_{j}-∑_{j=1}^{p}α_{j}+1 ρ=∏_{j=1}^{p}α_{j}^{-α_{j}}⋅∏_{j=1}^{q}β_{j}^{β_{j}} ϑ=½(q-p)+∑_{j=1}^{p}a_{j}-∑_{j=1}^{q}b_{j} If κ>0 the series has an infinite radius of convergence and _p\Psi_q(ζ) is an entire function. Series is uniformly and absolutely convergent for all finite ζ . If κ<0 the sum is divergent for all nonzero values of ζ whereas for κ=0 the function series has a finite radius of convergence ρ. Convergence on the boundary |ζ|=ρ depends on parameter ϑ converging absolutely if ℜ(ϑ)<-½.
For |arg(-ζ)|<π-ε, the Mellin-Barnes Integral _{p}Ψ_{q}(ζ)=\frac{1}{2πi}\int_{L}\Gamma(s)⋅\frac{\prod_{j=1}^{p}\Gamma(a_{j}-α_{j}s)}{\prod_{k=1}^{q}\Gamma(b_{k}-\beta_{k}s)}(-ζ)^{-s}ds defines a wider representation of Wright function. L is a contour separating the poles of \Gamma(s) to the left from those of \Gamma(a_{j}-α_{j}s) to the right. For contour L going from -i\infty up to +i\infty (possibly non-parallel to the vertical axis) this integral provides an analytical continuation of _{p}Ψ_{q}(ζ) in ζ∈ℂ\backslash [ρ,∞) when κ=0.
This function is a special case of FoxH function, (See Wiki's or Wolfram's sites)
_p\Psi_q(A;B;ζ)=H_{1+q,p}^{p,1}((1,1),B;A;-ζ^{-1}) For this particular case A=[(1,1),(1/a,2/a)] and B=[(1,2)]. Thus G_a function is G_a(x)=\frac{_2\Psi_1([(1,1),(1/a,2/a)];[(1,2)];x)}{\Gamma(1/a)} G_a(x)=\frac{H_{2,2}^{2,1}([[(1,1)],[(1,2)]];[[(1,1),(1/a,2/a)],[\cdot]];-x^{-1})}{\Gamma(1/a)} Note that κ=2(1-1/a) and series converges for a>1 to an entire function. You can set this expression using FoxH function in Wolfram's Mathematica v13.0 in symbolic mode to see if there are some explicit formulae for other values of parameter a. I suggest try with a\in \mathbb{Q} where a>1