# Compute the two-fold partial integral, where the three-fold full integral is known

I have the following trivariate ($$\rho_{11}, \rho_{22}, \mu$$) function $$$$4 \mu ^{3 \beta +1} \rho_{11}^{3 \beta +1} \left(-\rho_{11}-\rho_{22}+1\right){}^{3 \beta +1} \rho_{22}^{3 \beta +1} \left(\mu ^2 \rho_{22}+\rho_{11}\right){}^{-3 \beta -2},$$$$ the (three-fold) integral (for $$\beta$$ nonnegative integer) of which over $$\mu \in [0,1]$$, $$\rho_{11} \in [0,1]$$,$$\rho_{22} \in [0,1- \rho_{11}]$$ is $$$$\frac{\Gamma \left(\frac{3 \beta }{2}+1\right)^4}{\Gamma (6 \beta +4)}.$$$$ I would like to know the (two-fold) integral, say $$f(\mu,\beta)$$, over $$\rho_{11} \in [0,1]$$,$$\rho_{22} \in [0,1- \rho_{11}]$$. (The three-fold integral result is based on the application of the Mathematica FindSequenceFunction command, and does not constitute a formal proof--as remarks of C. Dunkl lead me to state.)

A reference for the (quantum-information-theoretic) background of this problem is sec. II.A.1 of Slater - Extended Studies of Separability Functions and Probabilities and the Relevance of Dyson Indices, in particular eq. (8) there. (For later related work, see secs. III,IV of Slater - Master Lovas-Andai and Equivalent Formulas Verifying the 8/33 Two-Qubit Hilbert-Schmidt Separability Probability and Companion Rational-Valued Conjectures, where I try to relate Lovas-Andai and Slater "separability functions".)

Mathematica computes the integral for specific nonnegative integer $$\beta$$, but apparently not for general $$\beta$$. (I can compute for $$\beta=1,2,\ldots$$ and then try to employ the Mathematica command FindSequenceFunction to obtain the general [two-fold] rule--which is, I think, how I got the three-fold result.)

For example, for $$\beta=1$$, we have $$$$\frac{\mu ^4 \left(12 \left(\mu^8+16\mu ^6+36 \mu ^4+16\mu ^2+1\right) \log (\mu )-5 \left(5 \mu ^8+32 \mu ^6-32 \mu ^2-5\right)\right)}{945 \left(\mu ^2-1\right)^9}$$$$ and for $$\beta=2$$, $$$$\frac{\mu ^7 \left(140 \left(\mu^2+1\right)\left(\mu ^{12}+48 \mu ^{10}+393 \mu ^8+832 \mu ^6+393 \mu ^4+48 \mu ^2+1\right)\log(\mu ) -(\mu^2-1) \left(363 \mu ^{12}+10310 \mu ^{10}+58673 \mu^8+101548 \mu ^6+58673 \mu ^4+10310 \mu ^2+363\right)\right) }{900900 \left(\mu ^2-1\right)^{15}}.$$$$

• I would ask for $$\int_{p=0}^1\int_{q=0}^{1-p}(\mu pq(1-p-q))^b(\mu^2q+p)^{-b-1}dp\ dq.$$ The messy variable names from the original are valuable context, but can be left out of the integral for clarity. – Matt F. Feb 11 '19 at 15:57
• It seems that the @MattF.'s integral is a linear polynomial in $\log\mu$ (I did not attempt to prove this), whose constant term has denominator $(\mu^2-1)^{2b}$ and whose linear term has denominator $(\mu^2-1)^{2b+1}$. – Martin Rubey Feb 11 '19 at 18:28
• The linear polynomials (in log $\mu$) that Martin Rubey observes appear to be multiples of $\mu ^{3 \beta +1} \left(\mu ^2-1\right)^{-6 \beta -3}$. – Paul B. Slater Feb 11 '19 at 19:28
• and moreover, the constant term in the linear term is the reciprocal of the $n$-th Catalan number, times $4(4b^2-1)$ – Martin Rubey Feb 11 '19 at 19:32
• In fact, the linear term (in @MattF.'s integral) is apparently $\frac{1}{(m^2-1)^{2b+1}}\frac{1}{4(4b^2-1)}\frac{b}{\binom{2(b-1)}{b-1}}\sum_k\binom{b}{k}^2m^{2k}$. – Martin Rubey Feb 11 '19 at 19:45

Using FriCAS it is actually not hard to guess the complete solution. However, this is not a proof, one would have to show

• that the integral is indeed a linear polynomial of $$\log\mu$$, (@Nemo?)
• that both the linear and the constant term of this polynomial must satisfy a second order recurrence with polynomial coefficient (@Nemo?),
• and that the degrees of the coefficients in the recurrence can be bounded.

It came as a surprise, but it probably shouldn't, that the recurrences for both the linear and he constant term are the same, they differ only in the initial conditions.

Anyway, here is code and solution. Be warned that I took out the $$\mu$$, and a factor of $$4$$ in the linear term.


f := (p*q*(1-p-q))^b*(m^2*q+p)^(-b-1)
)se fu ca all fint
fint(b:NNI):EXPR INT == integrate(integrate(eval(f, 'b=b), q=0..1-p, "noPole")::EXPR INT, p=0..1, "noPole")::EXPR INT
logm := first kernels fint(1)::EXPR INT
fintUP b == eval(fint(b), logm = lm)::UP(lm, FRAC POLY INT)

r1 := guessPRec([coefficient(fintUP b, 1) for b in 0..15], indexName=='b)
r2 := guessPRec([coefficient(fintUP b, 0) for b in 0..15], indexName=='b)

)expose RECOP
fun := operator 'f
r1DMP := eval(getEq(r1.1), [fun b = 'fb0, fun(b+1) = 'fb1, fun(b+2) = 'fb2])::DMP(['fb,'fb1,'fb2], POLY INT);
r2DMP := eval(getEq(r2.1), [fun b = 'fb0, fun(b+1) = 'fb1, fun(b+2) = 'fb2])::DMP(['fb,'fb1,'fb2], POLY INT);

factor coefficient(r1DMP, 'fb0, 1)
factor coefficient(r1DMP, 'fb1, 1)
factor coefficient(r1DMP, 'fb2, 1)


The recurrence for both the constant and the linear term is: $$-(b+1)^2 f(b) + 2 (2b+3)^2(\mu^2 +1) f(b+1) - 4 (2b+3)(2b+5)(\mu^2-1)^2 f(b+2) = 0$$

The initial conditions for the linear term are: $$f(0)=\frac{1}{2(\mu^2-1)},\quad f(1)=\frac{\mu^2+1}{12(\mu^2-1)^3}$$ The initial conditions for the constant term are: $$f(0)=0,\quad f(1)=\frac{1}{3(\mu^2-1)^2}$$

• That's great, Martin! I'll have to cogitate on this further in the A. M. First I heard of FriCAS. Can we get an explicit expression for the constant term akin to that I have for the linear (in log(\mu)) term? (I was getting pretty confident that the two terms had "a lot in common".) Also, you've been apparently using $m$ and $\mu$ for the same variable. Once things are clearly settled, I'll have to think about how to employ the results in the quantum-information-theoretic context that I referenced in my initial question. – Paul B. Slater Feb 14 '19 at 4:59
• So, Martin, I presume there is a summation formula for the constant term, parallel to that you gave for the linear term. Mathematica "automatically" converted that formula into the hypergeometric-based one. So, a parallel formula should give me the natural companion v(b,\mu) to my w(b,\mu) expression. – Paul B. Slater Feb 14 '19 at 13:13

The question as originally posed was to compute the integral of the trivariate function $$$$4 \mu ^{3 \beta +1} \rho_{11}^{3 \beta +1} \left(-\rho_{11}-\rho_{22}+1\right){}^{3 \beta +1} \rho_{22}^{3 \beta +1} \left(\mu ^2 \rho_{22}+\rho_{11}\right){}^{-3 \beta -2}.$$$$ Matt F. immediately suggested “streamlining” it to the form (transforming $$\beta$$ into $$\frac{b-1}{3}$$) $$$$\int_{p=0}^1 \int_{q=0}^{1-p} (\mu p q (1-p-q))^b (\mu^2 q +p)^{-b-1} dq dp,$$$$ which has been employed here since.

Subsequent computations, then, revealed that the two-fold integrals for $$b=1,2,\ldots$$ took the form $$$$v(b,\mu) + w(b,\mu) \log(\mu).$$$$ In a comment to the question, Rubey was able to express the coefficient $$w(b,\mu)$$ of the $$\log{\mu}$$ term as $$$$\frac{b \left(\mu ^2-1\right)^{-2 b-1} \sum _{k=0}^b \mu ^{2 k} \binom{b}{k}^2}{4 \left(4 b^2-1\right) \binom{2 (b-1)}{b-1}}$$$$ (note, in particular, the summation of a squared binomial). This required an additional factor of $$4 \mu^b$$, as well as a notational correction of $$m$$ to $$\mu$$ (as Rubey agreed in a subsequent comment).

Mathematica, interestingly, converted the Rubey expression (performing the indicated summation) to the hypergeometric-based formula, $$$$w(b,\mu)=\frac{\sqrt{\pi } 4^{-b} \mu^b \left(\mu^2-1\right)^{-2 b-1} \Gamma (b+1) \, _2F_1\left(-b,-b;1;\mu^2\right)}{\Gamma \left(b+\frac{3}{2}\right)}.$$$$

In his later answer to the question, Rubey showed that $$w(b,\mu)$$ and $$v(b,\mu)$$ were both given by the same second order recurrence, but with different initial conditions. Making use of these results we have found that $$$$v(b,\mu)= -\frac{\sqrt{\pi } 4^{-b} b! \mu ^b \left(\mu ^2-1\right)^{-2 b} \sum _{k=1}^b \mu ^{2 k-2} \sum _{i=0}^{k-1} \binom{b}{i}^2 (\psi ^{(0)}(b-i+1)-\psi ^{(0)}(i+1))}{\Gamma \left(b+\frac{3}{2}\right)},$$$$ where $$\psi$$ denotes the polygamma function. (The process of deriving this formula, in part, was the subject of my posting https://math.stackexchange.com/questions/3115582/do-these-polynomials-with-harmonic-number-related-coefficients-lie-in-some-parti .)

The $$v(b,\mu)$$ and $$w(b,\mu)$$ formulas have been explicitly found to yield the two-fold integration results for $$b=1,2,\ldots,11$$.

Most interestingly now, we observe that unlike $$w(b,\mu)$$, this formula for $$v(b,\mu)$$ does not appear readily to have an equivalent hypergeometric expression (despite what seem like simpler initial conditions for $$v(b,\mu)$$ in the second order recurrences). Both formulas involve the summation of a squared binomial, but the summand now has an added polygamma-based factor.

This phenomenon (failure to convert the summation) would appear to be a complicating factor in the underlying quantum-information-theoretic question we have been hoping to address. This involves finding “separability” functions $$f(b,\mu)$$ which when multiplied by the normalization (by $$\frac{\Gamma (b)^4}{\Gamma (2 b+2)}$$) of $$v(b,\mu) +w(b,\mu) \log{\mu}$$, and integrated over $$\mu \in [0,1]$$ would yield certain target (separability probability) values, given by $$$$1-\sqrt{\pi } 2^{-\frac{3 b}{2}-1} \Gamma \left(\frac{5 b}{12}+\frac{47}{24}\right) \Gamma \left(\frac{b}{2}+1\right) \Gamma \left(\frac{5 b}{6}+\frac{7}{6}\right) \Gamma \left(\frac{2 (b+2)}{3}\right) \, _6\tilde{F}_5\left(1,\frac{1}{6} (2 b+7),\frac{1}{12} (5 b+7),\frac{1}{12} (5 b+13),\frac{1}{24} (10 b+47),\frac{b+2}{2};\frac{b+11}{6},\frac{1}{24} (10 b+23),\frac{5 b}{12}+\frac{4}{3},\frac{1}{12} (5 b+22),\frac{2 (b+2)}{3};1\right)$$$$ (where the regularized hypergeometric function is indicated). For $$b=4, 7,10$$, this formula gives $$\frac{29}{64}$$, $$\frac{8}{33}$$ and $$\frac{26}{323}$$, the two-rebit, two-qubit and two-"quaterbit" Hilbert-Schmidt separability probabilities, respectively. https://arxiv.org/abs/1701.01973 (Any advice pertaining to this matter would be appreciated.)

Assembling the two functions of interest in the question posed, we can assert that $$$$v(b,\mu) + w(b,\mu) \log(\mu)=$$$$ $$$$\frac{1}{\Gamma \left(b+\frac{3}{2}\right)} \sqrt{\pi } 4^{-b} \mu^b \left(\mu^2-1\right)^{-2 b-1} \Gamma (b+1) \left(\log (\mu) \sum _{k=0}^b \mu ^{2 k} \binom{b}{k}^2-\left(\mu^2-1\right) \sum _{k=1}^b \mu^{2 k-2} \sum _{i=0}^{k-1} \binom{b}{i}^2 (\psi ^{(0)}(b-i+1)-\psi ^{(0)}(i+1))\right),$$$$ where, as previously noted, $$$$\sum _{k=0}^b \mu ^{2 k} \binom{b}{k}^2=\, _2F_1\left(-b,-b;1;\mu ^2\right).$$$$ (Perhaps the "master equation" just given can be somewhat "streamlined" in appearance--for instance, having all summations start with index zero.)

This completes our earlier "partial/half" answer to the question, which lacked the new interesting (less tractable, perhaps) expression for $$v(b,\mu)$$. My apologies for this somewhat awkward process of answering the posed question. (I considered simply replacing/editing the original answer, but that already had several comments pertaining to it--which might be somewhat confusing for subsequent comments.)

• Congratulations! In the last formula, you probably want to replace $u$ with $\mu$ throughout. Possibly it makes sense to write out the ${}_2F_1$, so the similarity between the two terms becomes more apparent. Finally, I assume that this is still a guess, or do you have a proof? @Nemo's method should provide a proof, though. – Martin Rubey Feb 21 '19 at 8:20
• Thanks again, Martin! Certainly, not a proof! Just observations of patterns for b=1,2.. and use of Mathematica FindSequenceFunction command (my frequent tool). The details of my reasoning process are in math.stackexchange.com/questions/3115582/… Perhaps the harmonic number aspect here should be stressed somewhat more. Can one employ the recurrences themselves directly in integrations over $\mu$ rather than having to find "closed"/hypergeometric forms for functions of $b$? How would the Nemo method work? – Paul B. Slater Feb 21 '19 at 15:42
• Regarding the remark of Martin Rubey that "@Nemo's method should provide a proof", I attempted unsuccessfully to do so in mathoverflow.net/questions/324114/… – Paul B. Slater Feb 26 '19 at 5:15

To begin, a partial/half answer based on the summation formula of Martin Rubey, given that the two-fold integration result in the Matt F. streamlined reformulation takes the form $$$$v(b,\mu) + w(b,\mu) \log(\mu) = \int_{p=0}^1 \int_{q=0}^{1-p} (\mu p q (1-p-q))^b (\mu^2 q +p)^{-b-1} dq dp,$$$$ is that (in line with the "hypergeometric" comment of Nemo above) $$$$w(b,\mu)=\frac{\sqrt{\pi } 4^{-b} \mu^b \left(\mu^2-1\right)^{-2 b-1} \Gamma (b+1) \, _2F_1\left(-b,-b;1;\mu^2\right)}{\Gamma \left(b+\frac{3}{2}\right)}.$$$$ The first term in the commented Rubey summation formula, $$\frac{b \left(\mu ^2-1\right)^{-2 b-1} \, _2F_1\left(-b,-b;1;\mu ^2\right)}{4 \left(4 b^2-1\right) \binom{2 (b-1)}{b-1}}$$, for $$w(b,\mu)$$ above (with $$m$$ corrected to $$\mu$$) required an additional factor of $$4 \mu^b$$ (as Ruhey noted in a subsequent comment). (Mathematica insists, it seems, on converting his summation to hypergeometric form. It would be of interest to know how Ruhey arrived at his formula--and if it might be of potential use in obtaining $$v(b,\mu)$$, as well.)

We have, additionally, trying to address the other "half" of the problem, found that $$$$v(b,1)=\frac{\pi 4^{-2 b-1} \Gamma (b+1)^2}{\Gamma \left(b+\frac{3}{2}\right)^2}.$$$$

Further, for $$v(b,\frac{1}{2})$$, Mathematica has given DifferenceRoot[Function[{[FormalY],[FormalN]},{(1+[FormalN])^2 [FormalY][[FormalN]]-5 (3+2 [FormalN])^2 [FormalY][1+[FormalN]]+(135+144 [FormalN]+36 [FormalN]^2) [FormalY][2+[FormalN]]==0,[FormalY][1]==-(8/27),[FormalY][2]==-(8/81)}]][b]

The TeX version of this that Mathematica gives is not accepted here because the "Argument to unicode must be a number". (This appears to be a problem with the TeXForm command of Mathematica https://mathematica.stackexchange.com/questions/191414/nonacceptance-by-stackexchange-site-of-mathematica-texform-employing-unicode .)

Along similar lines, for $$v(b,1/3)$$, we have DifferenceRoot[ Function[{[FormalY], [FormalN]}, {9 (1 + [FormalN])^2 \ [FormalY][[FormalN]] - 60 (3 + 2 [FormalN])^2 [FormalY][ 1 + [FormalN]] + (3840 + 4096 [FormalN] + 1024 [FormalN]^2) [FormalY][2 + [FormalN]] == 0, [FormalY][1] == -(9/64), [FormalY][2] == -(81/4096)}]][b]

So, "DifferenceRoot" formulas appear to be the case for other specific values of $$v(b,\mu)$$ with $$\mu \neq 1$$.

Perhaps, we will attempt ("manually") to present the DifferenceRoot formulas above in more standard, proper TeX form (per comment of Somos in https://mathematica.stackexchange.com/questions/191414/nonacceptance-by-stackexchange-site-of-mathematica-texform-employing-unicode ).

It may be helpful to note that the inner integral of the Matt F. reformulation can be performed, yielding $$$$\int_{q=0}^{1-p} (\mu p q (1-p-q))^b (\mu^2 q +p)^{-b-1} dq =-\frac{(p-1) \Gamma (b+1)^2 \left((p-1)^2 \mu\right)^b \, _2\tilde{F}_1\left(b+1,b+1;2 (b+1);\frac{(p-1) \mu^2}{p}\right)}{p},$$$$ where the regularized hypergeometric function is indicated.

So, we would like the counterpart for $$v(b,\mu)$$ of the Rubey formula for $$w(b,\mu)$$, that is, $$$$4 u^b \left(u^2-1\right)^{-2 b-1} \frac{1}{4 \left(4 b^2-1\right)} \frac{b} {\binom{2 (b-1)}{b-1}} \Sigma_{k=0}^b u^{2 k} \binom{b}{k}^2,$$$$ which he apparently obtained using the general purpose computer algebra system, FriCAS.

A new, interesting observation is that if we multiply $$v(b,\mu)$$ by $$$$(-1+\mu^2)^{2 b+1} \mu^{-b} P,$$$$ where, with $$\gamma$$ being Euler's constant, $$\psi ^{(0)}(b+1)$$ being the zero-th derivative of the digamma function, and $$H_b$$, the $$b$$-th harmonic number, $$$$P= \frac{\sqrt{\pi } 2^{-2 b} b! (\psi ^{(0)}(b+1)+\gamma )}{\left(b+\frac{1}{2}\right)!} = \frac{\sqrt{\pi } 4^{-b} b! H_b}{\Gamma \left(b+\frac{3}{2}\right)},$$$$ we obtain ("asymmetric"-type) polynomials in $$\mu$$ of degree $$2 b$$, all with constant terms equal to 1 and highest terms equal to $$-u^{2 b}$$. For example, for $$b=1$$, we have $$$$1-\mu^2$$$$ for $$b=2$$, $$$$1-\mu^4$$$$ for $$b=3$$, $$$$-\mu ^6-\frac{27 \mu ^4}{11}+\frac{27 \mu ^2}{11}+1$$$$ for $$b=4$$, $$$$-\mu ^8-\frac{32 \mu ^6}{5}+\frac{32 \mu ^2}{5}+1,$$$$ for $$b=5$$, $$$$-\mu ^{10}-\frac{1625 \mu ^8}{137}-\frac{2000 \mu ^6}{137}+\frac{2000 \mu ^4}{137}+\frac{1625 \mu ^2}{137}+1, \ldots$$$$ So, one needs to find the governing rule for these polynomials, which would complete the formula for $$v(b,\mu)$$.

• I am quite sure that I do get the factor $(\mu^2-1)^{2b+1}$, could you please recheck? For example, for $b=2$, I obtain $$\frac{\mu^4+4\mu^2+1}{15(\mu^2-1)^5}$$. – Martin Rubey Feb 13 '19 at 6:59
• Thanks, Martin! I did some quick checking this morning--before off to campus--and I can't, at the moment, get anything to agree. Seems that I may need a factor of u^b/(u^2-1)^(2 b+1) in my formula. Embarrassing! – Paul B. Slater Feb 13 '19 at 14:44
• Oh yes, I took out the $\mu^b$ from the integral! – Martin Rubey Feb 13 '19 at 15:29
• So, the revised answer accounts for these last matters. – Paul B. Slater Feb 13 '19 at 22:37
• Matt F. I think the result (according to Mathematica) should be $\, _2F_1\left(-b,-b;1;u^2\right)$. $k$ is the index of summation and should not appear in the answer. – Paul B. Slater Feb 14 '19 at 21:48

The most concise answer to the problem (in the Matt F. streamlined form) $$$$\int_{p=0}^1 \int_{q=0}^{1-p} (\mu p q (1-p-q))^b (\mu^2 q +p)^{-b-1} dq dp,$$$$ so far given is $$$$\frac{\pi (-1)^{2 b} 4^{-2 b-1} \mu^{-b-2} \Gamma (b+1)^2 \, _2F_1\left(b+1,b+1;2 (b+1);1-\frac{1}{\mu^2}\right)}{\Gamma \left(b+\frac{3}{2}\right)^2}.$$$$ Its validity can be seen by, first, performing the inner integration of $$q \in [0,1-p]$$, obtaining thereby, $$$$\int_{p=0}^1-\frac{(p-1)^{2 b+1} \mu^b \Gamma (b+1)^2 \, _2F_1\left(b+1,b+1;2 (b+1);\frac{(p-1) \mu^2}{p}\right)}{p \Gamma (2 (b+1))} dp.$$$$ Then, the hypergeometric argument $$\frac{(p-1) \mu^2}{p}$$ is transformed to $$v$$, by performing the change-of-variables $$p\to \frac{\mu^2}{\mu^2-v}$$. The one-fold integration problem, which can be carried out, then, becomes
$$$$\int_{v=-\infty}^0 \mu ^b \left(-v^{2 b+1}\right) \Gamma (b+1)^2 \, _2\tilde{F}_1(b+1,b+1;2 b+2;v) \left(\mu ^2-v\right)^{-2 (b+1)} dv,$$$$ yielding us the indicated answer given above here to the two-fold integration (in the Matt F. form).

However, to perform the third step of the three-fold integration, that is, $$$$\int_{\mu=0}^1 \frac{\pi (-1)^{2 b} 4^{-2 b-1} \mu^{-b-2} \Gamma (b+1)^2 \, _2F_1\left(b+1,b+1;2 (b+1);1-\frac{1}{\mu^2}\right)}{\Gamma \left(b+\frac{3}{2}\right)^2} d \mu,$$$$ we first perform--as pointed out by Charles Dunkl--the quadratic transformation 15.8.13 given in https://dlmf.nist.gov/15.8#iii . That is $$$$\, _2F_1(a,b;2 b;z) = (1-\frac{1}{2} z)^{-a} \, _2F_1(\frac{1}{2} a,\frac{1}{2} a+\frac{1}{2};b+\frac{1}{2} ;(\frac{z}{2-z})^2)$$$$ on the integrand and make the change-of-variables $$\mu\to \sqrt{\frac{2}{\sqrt{Y}+1}-1}$$ and now integrate over $$Y \in [1,0]$$, that is, $$$$\int_{Y=1}^0 \frac{\pi 2^{-4 b-3} e^{2 i \pi b} \left(\frac{2}{\sqrt{Y}+1}-1\right)^{\frac{1-b}{2}} \left(\frac{1}{1-\sqrt{Y}}\right)^{-b} \Gamma (b+1)^2 \, _2F_1\left(\frac{b+1}{2},\frac{b+2}{2};b+\frac{3}{2};Y\right)}{\left(Y-\sqrt{Y}\right) \Gamma \left(b+\frac{3}{2}\right)^2} dY.$$$$ The result of this integration is $$$$\frac{e^{2 i \pi b} \Gamma \left(\frac{b+1}{2}\right)^4}{4 \Gamma (2 b+2)}.$$$$ For integral $$b$$, the phase factor $$e^{2 i \pi b}$$ is, of course, 1, and for half-integral $$b$$, it is -1.

Transforming back to the original (Dyson-index-like parameter) $$\beta$$, using $$b \to 1+3 \beta$$, we obtain $$$$\frac{e^{6 i \pi \beta } \Gamma \left(\frac{3 \beta }{2}+1\right)^4}{4 \Gamma (6 \beta +4)}.$$$$

This agrees (up to the omitted--in the Matt F. parameterization--factor 4) with the original assertion of $$$$\frac{\Gamma \left(\frac{3 \beta }{2}+1\right)^4}{\Gamma (6 \beta +4)},$$$$ for integral $$\beta$$.

Here is a link--through the wolframcloud--to a Mathematica notebook presenting the main later stages of our argument.

https://www.wolframcloud.com/objects/e8bd5d73-06a4-4798-8485-8e2597f660eb