SE overestimation for extreme p-values due to incorrect log-p handling in R::qnorm call (regression in C++ rewrite)

[DDochtermann]

There is a bug in src/SAIGE_test.cpp in the SE back-calculation after SPA is applied, affecting markers where the score-test p-value underflows to 0 (i.e. large N, common variants with strong effects). This appears to be a regression of the same class of bug documented in the v1.3.0 release notes ("fix the bug that causes overestimation of SE for markers when SPA is applied... introduced since v1.1.6"), which was fixed on the R side but re-introduced in the C++ implementation.

src/SAIGE_test.cpp, line 662, inside SAIGEClass::getMarkerPval:
t_qval = R::qnorm(t_SPApval/2, 0, 1, false, ispvallog);

ispvallog is set to true in the score test functions when the chi-squared p-value underflows to exactly 0 (i.e. boost::math::cdf returns 0 for extreme statistics). In that case, t_SPApval holds a log-scale p-value (a large negative number). The call to R::qnorm correctly passes ispvallog as the log_p argument, but the first argument is always computed as t_SPApval/2, which is only correct when t_SPApval is on the linear scale. When ispvallog=true, the correct log-scale operation is t_SPApval - log(2), not /2.

As a result, R::qnorm receives a nonsensical value as its probability argument, t_qval is wrong, and t_seBeta = fabs(t_Beta)/t_qval is overestimated for all markers that triggered the log-p path.

Likely fixed with
t_qval = R::qnorm(ispvallog ? t_SPApval - std::log(2.0) : t_SPApval/2, 0, 1, false, ispvallog);

Any variant where the score-test chi-sq statistic is large enough that boost::math::cdf(complement(chisq_dist, stat)) returns exactly 0 triggers the issue. In large biobank cohorts (e.g. N > 400k) this occurs routinely for well-powered common variants. The affected rows are identifiable in output as those where p.value and p.value.NA differ and p.value is near the double-precision floor (~1e-306 range).

For affected markers, SE is overestimated, and 2*pnorm(-abs(BETA/SE)) recovers p.value.NA (the pre-SPA score test p) rather than p.value (the SPA p). Markers where the score-test p does not underflow are unaffected.

The v1.3.0 changelog notes that the SE overestimation bug from v1.1.6 was fixed, with the documented correction being SE = abs(BETA/qnorm(p.value/2)). That fix addressed the R-layer code. The same error is present in the C++ reimplementation of the same logic.