code: def estimate_chi2mixture(self, lrt):
          
  
          """
          step 1: estimate the probability of being in component one
          """
          self.mixture = 1-(lrt<=self.tol).mean()
          n_false      = SP.sum(lrt>self.tol)
  
          """
          step 2: only use the largest qmax fraction of test statistics to estimate the
                  remaining parameters
          """
          n_fitting   = SP.ceil(self.qmax * n_false)
          lrt_sorted  = -SP.sort(-lrt)[:n_fitting]
          q           = SP.linspace(0, 1,n_false)[1:n_fitting+1]
          log_q       = SP.log10(q)
          
          """
          step 3: fitting scale and dof by minimizing the squared error of the log10 p-values
                  with their theorietical values [uniform distribution]
          """
          MSE_opt = SP.inf
          MSE     = SP.zeros((self.n_intervals,self.n_intervals))
  
          
          for i,scale in enumerate(SP.linspace(self.scale_min,self.scale_max,self.n_intervals)):
              for j,dof in enumerate(SP.linspace(self.dof_min,self.dof_max,self.n_intervals)):
                  p     = STATS.chi2.sf(lrt_sorted/scale,dof)
                  log_p = SP.log10(p)
                  MSE[i,j]   = SP.mean((log_q - log_p)**2)
                  if MSE[i,j] < MSE_opt:
                      MSE_opt    = MSE[i,j]
                      self.scale = scale
                      self.dof   = dof
docs: estimates the parameters of a mixture of a chi-squared random variable of degree
          0 and a scaled chi-squared random variable of degree d
          
          (1-mixture)*chi2(0) + (mixture)*scale*chi2(dof),
  
          where
          scale   is the scaling parameter for the scales chi-square distribution
          dof     are the degrees of freedom of the second component
          mixture is the probability of beeing in the first component
  
          input:
          lrt         [Ntests] vector of test statistics