qqplot and example

extended_line.png

qqplot.ml

let _draw_extended_line ?(h=_default_handle) x0 y0 x1 y1 l r u d =
    (*Specify two points, draw a line between them, and extend on both ends with dash line*)
    let x0, y0, x1, y1 = if x0 > x1 then x1, y1, x0, y0 else x0, y0, x1, y1 in
    let yl = if x0 = x1 then u else y0 -. (y1 -. y0) /. (x1 -. x0) *. (x0 -. l) in
    let yr = if x0 = x1 then d else y1 +. (y1 -. y0) /. (x1 -. x0) *. (r -. x1) in
    let xl = if x0 = x1 then x0 else l in
    let xr = if x0 = x1 then x0 else r in
    let _ = draw_line ~h ~line_style:1 x0 y0 x1 y1 in
    let _ = draw_line ~h ~line_style:3 xl yl x0 y0 in
    draw_line ~h ~line_style:3 x1 y1 xr yr
;;

let qqplot ?(h=_default_handle) ?(color=(-1,-1,-1)) ?(marker_size=4.)
  ?(pd=(fun i -> Owl_stats.Cdf.gaussian_Pinv i 1.)) ?y x =
  (*Note: data in x are actually plotted along y-axis*)
  (*Input: x and y should be vectors (to be extend to matrices of same col number)
    with more than one element; their lengths should be equal (a condition to be removed) *)
  let x = Owl_dense_matrix.D.to_array x |> Owl_stats.sort ~inc:true in
  let y = match y with
    | Some arr -> Owl_dense_matrix.D.to_array arr |> Owl_stats.sort ~inc:true
    | None     -> let n = Array.length x in
                    let qth = Owl_dense_matrix.D.linspace ((1. -. 0.5) /. float_of_int n) (( float_of_int n-. 0.5) /. float_of_int n) n in
                    let q = Owl_dense_matrix.D.map pd qth in
                    Owl_dense_matrix.D.to_array q
  in
  (*A line joining the first and third quartiles of each distribution*)
  let p1y, p1x = (Owl_stats.first_quartile y, Owl_stats.first_quartile x) in
  let p3y, p3x = (Owl_stats.third_quartile y, Owl_stats.third_quartile x) in
  let l, r = Owl_stats.minmax y in
  let u, d = Owl_stats.minmax x in
  let x = Owl_dense_matrix.D.of_array x 1 (Array.length x) in
  let y = Owl_dense_matrix.D.of_array y 1 (Array.length y) in
  let _ = scatter ~h ~color:color ~marker:"#[0x002b]" ~marker_size:marker_size y x in
  _draw_extended_line ~h p1y p1x p3y p3x l r u d
;;
  
(*Example*)
let x = Mat.(gaussian 100 1 *$ 10.) in 
let y1 = Mat.gaussian 100 1 in
let h = create ~m:2 ~n:2 "" in
let _ = set_background_color h 28 28 28 in
let _ = subplot h 0 0 in
let _ = set_title h "Gaussian vs. Gaussian " in 
let _ = set_ylabel h "Quantiles of Input Sample" in 
let _ = set_xlabel h "Normal Distribution Quantiles" in 
let _ = qqplot ~h x ~y:y1  in 
let _ = subplot h 0 1 in
let _ = set_title h "Gaussian vs. Default Distribution" in 
let _ = set_ylabel h "Quantiles of Input Sample" in 
let _ = set_xlabel h "Normal Distribution Quantiles" in 
let _ = qqplot ~h x ~color:(0,128,255) in
let _ = subplot h 1 0 in
let _ = set_title h "Gaussian vs. Rayleigh" in 
let _ = set_ylabel h "Quantiles of Input Sample" in 
let _ = set_xlabel h "Rayleigh Distribution Quantiles" in 
let _ = qqplot ~h x ~pd:(fun p -> Stats.Cdf.rayleigh_Pinv p 0.5)  in
let _ = subplot h 1 1 in
let _ = set_title h "Gaussian vs. Chi-Sqare" in 
let _ = set_ylabel h "Quantiles of Input Sample" in 
let _ = set_xlabel h "Chi-Sqare Distribution (k=6) Quantiles" in 
let _ = qqplot ~h x ~pd:(fun p ->  Stats.Cdf.chisq_Pinv p 6.) ~marker_size:3. in
output h;;

(*Example: Extended Line*)
let h = create "" in
let _ = set_background_color h 28 28 28 in
let _ = _draw_extended_line ~h (-1.4) (-1.4) 1.4 1.4 (-4.) 4. 4. (-4.) in
let _ = _draw_extended_line ~h (-2.) 0. 2. 0. (-4.) 4. 4. (-4.) in 
let _ = _draw_extended_line ~h 0. 2. 0. (-2.) (-4.) 4. 4. (-4.) in 
let _ = _draw_extended_line ~h 1.4 (-1.4) (-1.4) 1.4 (-4.) 4. 4. (-4.) in 
output h;;

qqplot.pdf

really cool!