v1.7.0 — pprodnormal3() correctness fix
Pre-releaseFixes a correctness bug in pprodnormal3() that returned silently wrong probabilities on ill-conditioned covariance. If you use three-variable products (serial indirect effects a1 * a2 * b), install this now.
remotes::install_github("Data-Wise/RMediation", ref = "dev")Pre-release.
mainremains at 1.6.1, which is what CRAN serves. CRAN's update cadence window opens 2026-08-21; this release exists so the fix is available before then.
What was wrong
pprodnormal3() (and its alias p_prod3()) returned materially wrong values whenever the (X, Y) covariance block was ill-conditioned — with no warning, no error, and no diagnostic. Tightening tol did not help: the integrator converged to a stable but incorrect value.
cdf(), dist_quantile(), confint() and ci() invert this CDF by root-finding, so all of them inherited the error.
At rho = 0.999, across the support:
| q | 1.6.1 | 1.7.0 | reference |
|---|---|---|---|
| -2.0 | 0.0044691 (92% rel. error) | 0.05711903 | 0.05711903 |
| 0.0 | 0.3445185 | 0.49727060 | 0.49727060 |
| 0.5 | 0.6109671 (24%) | 0.80671502 | 0.80671502 |
| 2.0 | 0.6506248 | 0.92056384 | 0.92056384 |
Maximum error is now 5.9e-09, down from 2.7e-01.
A second, previously unreported failure is also fixed: at large standardized means (mean / sd beyond roughly 8) the old integrator returned exactly 0 where the true probability is around 0.46. That one was independent of conditioning.
What changed
The default integration method is now "gauss" — tensor-product Gauss-Legendre on a domain partitioned at the coordinate axes, with the node count escalated automatically until successive rules agree to tol.
- Results change. On well-conditioned input the new default agrees with the old to about 1e-7. On ill-conditioned input it differs substantially; that difference is the correction.
- Faster, by roughly 8-20x at the covariances tested.
method = "hcubature"still selects the previous integrator, for cross-checking. Not recommended otherwise.- New
nodesargument forces a fixed rule; newdiagnosticsargument exposes the convergence estimate.diagnosticsdefaults toFALSE, so the return value stays a bare numeric. pprodnormal3()now warns rather than returning a wrong answer silently when quadrature hits the node cap without meetingtol.
Verification
- Test suite: 339 passed, 0 failed, 0 warnings, 0 skipped
R CMD check --as-cran: 0 errors, 0 warnings, 0 notes- CI: 6/6 green across Ubuntu (release + devel), macOS, and Windows
Reference values are Gauss-Legendre results self-consistent to 1e-14 across n = 512 / 1024 / 2048, each independently within 3 SE of Monte Carlo at 1e7 draws.
Full detail: SPEC-27-prod3-integrator.md (design and evidence) and SPEC-28-review-triage.md (review triage).