machenko_pastur

[spinkney]

Also see https://www.wolframalpha.com/input/?i=MarchenkoPasturDistribution

https://github.com/JuliaStats/Distributions.jl/pull/1288/files

functions {
  marchenko_pastur_pdf(real x, real Q, real sigma) {
        real q = inv_sqrt(Q);
        real b = square(sigma * (1 + q)) # Largest eigenvalue
        real a = square(sigma * (1 - q)) # Smallest eigenvalue
        
        if (x > b || x < a) return 0;
         else return ( Q * sqrt((b - x) * (x - a)) ) / ( 2 * pi() * square(sigma) * x );
  }
  
    marchenko_pastur_lpdf(real x, real Q, real sigma) {
        real q = inv_sqrt(Q);
        real b = square(sigma * (1 + q)) # Largest eigenvalue
        real a = square(sigma * (1 - q)) # Smallest eigenvalue
        
        if (x > b || x < a) return negative_infinity();
         else return log(Q) + 0.5 * ( (b - x) * (x - a) ) - 2 * sigma - log(x); //  -log(2) - log(pi()) constants removed
  }

  marchenko_pastur_cdf(real x, real Q, real sigma) {
      real b = square(sigma * (1 + q)) # Largest eigenvalue
      real a = square(sigma * (1 - q)) # Smallest eigenvalue
      real R = (b - x) * (x - a);
      
      return sqrt(R) + 
             0.5 * (a + b) * asin( (2 * x - a - b) / (b - a) ) -
             ( a * b / sqrt(a * b) ) * asin( (a + b) * x - (2 * a * b) / (x * (b - a)) );
  }
}