Warrior_EA/System/BinomialStats.mqh

//+------------------------------------------------------------------+
//|                                                  BinomialStats.mqh |
//|                                                      AnimateDread  |
//| project. Free functions, no state, so the deploy gate, the two   |
//| edge floors, the barrier ladder and the detectability reports    |
//| all read the same formula instead of nine transcriptions of it.  |
//+------------------------------------------------------------------+
#ifndef WARRIOR_SYSTEM_BINOMIALSTATS_MQH
#define WARRIOR_SYSTEM_BINOMIALSTATS_MQH
#include <Math\Stat\Normal.mqh>
//--- Beta.mqh includes only Math.mqh, already pulled in by Normal.mqh above - no new transitive
//--- dependency. Guards only a<=0||b<=0, so it takes FRACTIONAL shape parameters natively, unlike
//--- MathCumulativeDistributionBinomial (Math\Stat\Binomial.mqh), which rejects non-integer n
//--- outright - and this file's effN is deflated for label overlap, essentially never an integer.
#include <Math\Stat\Beta.mqh>
//+------------------------------------------------------------------+
//| Upper-tail standard normal, Q(z) = P(Z >= z).                     |
//|                                                                   |
//| Wraps the library so the NaN policy lives in one place: an        |
//| unusable z reads as "not significant" rather than propagating a   |
//| NaN into a gate decision. tail=false asks for the UPPER tail, and |
//| the clamp keeps a -1e-17 round-off out of the Sidak power.        |
//+------------------------------------------------------------------+
double NormalUpperTailQ(const double z)
  {
   if(!MathIsValidNumber(z))
      return(1.0);
   int    err=0;
   double q=MathCumulativeDistributionNormal(z,0.0,1.0,false,false,err);
   if(err!=ERR_OK || !MathIsValidNumber(q))
      return(1.0);
   return(MathMax(0.0,MathMin(1.0,q)));
  }
//+------------------------------------------------------------------+
//| Sampling variance of a binomial rate p over n observations:       |
//| p(1-p)/n, in fraction^2.                                          |
//|                                                                   |
//| Returns 0 for a degenerate rate or an empty sample, which every   |
//| caller already treats as "no bar to clear". Callers that combine  |
//| symbols by inverse variance want this rather than the SE.         |
//+------------------------------------------------------------------+
double BinomialVar(const double p,const double n)
  {
   if(!MathIsValidNumber(p) || !MathIsValidNumber(n))
      return(0.0);
   if(n<=0.0 || p<=0.0 || p>=1.0)
      return(0.0);
   return(p*(1.0-p)/n);
  }
//+------------------------------------------------------------------+
//| Standard error of a binomial rate p over n observations, in       |
//| percentage points: 100 * sqrt(p(1-p)/n).                          |
//+------------------------------------------------------------------+
double BinomialSEPct(const double p,const double n)
  {
   return(100.0*MathSqrt(BinomialVar(p,n)));
  }
//+------------------------------------------------------------------+
//| EXACT one-sided binomial upper-tail p-value, P(Binomial(n,p0)     |
//| >= k), via the regularized-incomplete-beta identity               |
//| P(X>=k) = I_p0(k, n-k+1) = CDF of Beta(k, n-k+1) evaluated at p0.  |
//| Replaces the normal approximation (chance + sigmas*SE) the deploy |
//| gate used to compare precPct against directly - that approximation|
//| is anticonservative near the ~60-observation effN this gate       |
//| actually sees (no continuity correction, understates the tail),   |
//| by roughly 1pp of win rate at n~63, p0=0.5. k/n may be FRACTIONAL  |
//| (effN is deflated for label overlap) - this is then a principled  |
//| interpolation of the exact test, not a literal one, but it is the |
//| same interpolation the rest of this gate's arithmetic already     |
//| makes by treating effN as a real number everywhere else.          |
//+------------------------------------------------------------------+
double BinomialUpperTailP(const double k,const double n,const double p0)
  {
   if(!MathIsValidNumber(k) || !MathIsValidNumber(n) || !MathIsValidNumber(p0))
      return(1.0);                   // fail closed - a NaN must never read as "significant"
   if(n<=0.0 || k<=0.0)
      return(1.0);                   // no evidence at all
   if(k>=n)
      return(0.0);                   // every trial succeeded - always significant; guards b<=0 below
   if(p0<=0.0)
      return(0.0);                   // any success at all beats a zero-probability null
   if(p0>=1.0)
      return(1.0);                   // nothing beats a certain null
   //--- tail=TRUE, and the difference is the whole test. Beta.mqh documents this flag as "Flag to
   //--- calculate lower tail", so tail=false returns 1-I_p0 - the OPPOSITE of the identity above.
   //--- NormalUpperTailQ() above passes false because it genuinely wants Q(z)=1-CDF; this wants the
   //--- CDF itself, and copying that call's flag shipped a gate whose floor was 100.0% for every
   //--- input (bisection never sees a significant mid, so `hi` never leaves its 100.0 seed), which
   //--- BarUnreachable() then read as "no win rate can ever clear this" on every model and every
   //--- ensemble. Caught 2026-08-25 in the live log, not by the compiler and not by inspection.
   int err=0;
   double p=MathCumulativeDistributionBeta(p0,k,n-k+1.0,true,false,err);
   if(err!=ERR_OK || !MathIsValidNumber(p))
      return(1.0);                   // fail closed, same policy as NormalUpperTailQ
   return(MathMax(0.0,MathMin(1.0,p)));
  }
//+------------------------------------------------------------------+
//| EXACT analogue of "chance + sigmas*SE": the smallest observed     |
//| rate (a PERCENTAGE) whose BinomialUpperTailP() clears the         |
//| one-sided normal tail at `sigmas`, at this effN. Found by          |
//| bisection - the incomplete beta has no closed-form inverse in     |
//| this direction, and this runs once per era-end gate check, not    |
//| per tick, so 60 steps (~1e-16pp resolution) costs nothing.         |
//|                                                                    |
//| MONOTONE BY CONSTRUCTION: raising the candidate rate raises k,     |
//| which can only lower or hold P(X>=k) - so "smallest rate that      |
//| clears the target p-value" is a well-posed bisection on a          |
//| decreasing function.                                              |
//+------------------------------------------------------------------+
double ExactEdgeFloorPct(const double chancePct,const double effN,const double sigmas)
  {
   //--- effN<=0 -> the floor IS chancePct (zero evidence, zero SE), matching what the normal
   //--- approximation this replaces did (BinomialSEPct returns 0 for n<=0, so chance+sigmas*0 ==
   //--- chance) - NOT 100.0/"impossible". A genuinely unreachable bar is what the bisection below
   //--- converges to on its own when even a 100% observed rate fails to clear the target p-value.
   if(!MathIsValidNumber(effN) || effN<=0.0)
      return(MathMax(0.0,chancePct));
   if(!MathIsValidNumber(chancePct))
      return(100.0);                 // no reference rate to test against - fail toward "unreachable"
   double targetAlpha=NormalUpperTailQ(sigmas);
   double p0=MathMax(0.0,MathMin(1.0,chancePct/100.0));
   double lo=chancePct, hi=100.0;
   for(int i=0;i<60;i++)
     {
      double mid=0.5*(lo+hi);
      double k=mid/100.0*effN;
      if(BinomialUpperTailP(k,effN,p0)<=targetAlpha)
         hi=mid;   // already significant at mid - the floor is at or below it
      else
         lo=mid;   // not yet significant - the floor is above it
     }
   return(hi);
  }
//+------------------------------------------------------------------+
//| Independent observations needed to certify an edge of `edge` over |
//| a base rate p, at `sigmas` standard errors:                       |
//|   n = sigmas^2 * p(1-p) / edge^2                                  |
//|                                                                   |
//| BinomialSEPct solved for n. Answers "could this configuration     |
//| EVER prove an edge this size" - a property of the geometry, the   |
//| horizon and the window, which no amount of training moves.        |
//+------------------------------------------------------------------+
double BinomialCallsForEdge(const double p,const double edge,const double sigmas)
  {
   if(edge<=0.0 || p<=0.0 || p>=1.0)
      return(0.0);
   return(sigmas*sigmas*p*(1.0-p)/(edge*edge));
  }
//+------------------------------------------------------------------+
//| A rate shrunk toward a prior - the estimator, where everything   |
//| above is the test.                                               |
//+------------------------------------------------------------------+
double ShrunkRatePct(const double hits,const double n,const double priorPct,const double priorN)
  {
   bool havePrior=(priorN>0.0 && MathIsValidNumber(priorPct) && priorPct>=0.0);
   if(!MathIsValidNumber(hits) || !MathIsValidNumber(n) || n<=0.0)
      return(havePrior ? priorPct : 0.0);   // no evidence => the prior IS the estimate
   if(!havePrior)
      return(100.0*hits/n);
   return((hits+priorN*(priorPct/100.0))*100.0/(n+priorN));
  }
//+------------------------------------------------------------------+
//| Sidak family-wise p for the best of nTried candidates:            |
//| 1 - (1 - p_single)^N.                                             |
//|                                                                   |
//| The null of the MAXIMUM, not of a single draw. At the magnitudes  |
//| in play (p ~ 1e-4..1e-2, N ~ 10..1000) plain double precision is  |
//| ample - no need for the log1p/expm1 form MQL5 would not give us   |
//| anyway.                                                           |
//+------------------------------------------------------------------+
double SidakFamilyP(const double zObs,const int nTried)
  {
   return(1.0-MathPow(1.0-NormalUpperTailQ(zObs),(double)MathMax(nTried,1)));
  }
//+------------------------------------------------------------------+
//| Add-one-smoothed Monte-Carlo permutation p-value:                 |
//| (1 + atLeast) / (draws + 1), where atLeast is the count of null   |
//| draws at least as extreme as the observed statistic.              |
//|                                                                   |
//| Returns 1.0 for draws<=0 (an abandoned or empty null), matching   |
//| every call site's own guard - a truncated null is not a smaller   |
//| null, it is a wrong one.                                          |
//+------------------------------------------------------------------+
double PermutationPValue(const int atLeast,const int draws)
  {
   if(draws<=0)
      return(1.0);
   return((double)(1+atLeast)/(draws+1));
  }
//+------------------------------------------------------------------+
#endif