Even if you cannot find exactly equivalent gamma functions you ought to be able to translate gb
with SciPy's integration and root finding functions. The functions required can be obtained, e.g. (illustrating with some demo values)

For example

As you can see, the code constructed from the more basic functions produces the same answer, albeit more slowly.
Code
gamma[z_] := \!\(
\*SubsuperscriptBox[\(\[Integral]\), \(0\), \(\[Infinity]\)]\(
\*SuperscriptBox[\(t\), \(z - 1\)]
\*SuperscriptBox[\(E\), \(-t\)] \[DifferentialD]t\)\)
gamma[a_, z0_, z1_] := \!\(
\*SubsuperscriptBox[\(\[Integral]\), \(z0\), \(z1\)]\(
\*SuperscriptBox[\(t\), \(a - 1\)]
\*SuperscriptBox[\(E\), \(-t\)] \[DifferentialD]t\)\)
gammaregularized[a_, z1_] := gamma[a, 0, z1]/gamma[a]
cdf[\[Beta]_, p_] :=
Piecewise[{{gammaregularized[1/\[Beta], p/\[Beta]], p > 0}}]
p = 1.2;
FindRoot[cdf[\[Beta], p] - .95, {\[Beta], 1, If[p == 1, 1.1, p]}]
FindRoot[CDF[GammaDistribution[1/\[Beta], \[Beta]],
p] - .95, {\[Beta], 1, If[p == 1, 1.1, p]}]