violin plot

violin.rkt

#lang racket/base

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; provide

(provide violin)

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; require

(require math/statistics
         plot
         plot/utils
         racket/dict
         racket/function)

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; functions

(define (violin vals
                #:bandwidth [bandwidth (silverman vals)]
                #:x-min [x-min #f]
                #:x-max [x-max #f]
                #:y-min [y-min #f]
                #:y-max [y-max #f]
                #:color [color (interval-color)]
                #:style [style (interval-style)]
                #:line-color [line-color (interval-line1-color)]
                #:line-width [line-width (interval-line1-width)]
                #:line-style [line-style (interval-line1-style)]
                #:alpha [alpha (interval-alpha)]
                #:label [label #f])
  (define-values (f low high)
    (kde vals bandwidth))
  (define x-axis (const 0))
  (define x-min* (or x-min low))
  (define x-max* (or x-max high))
  (define settings
    `([#:y-min . ,y-min]
      [#:y-max . ,y-max]
      [#:color . ,color]
      [#:style . ,style]
      [#:line1-style . transparent]
      [#:line2-color . ,line-color]
      [#:line2-width . ,line-width]
      [#:line2-style . ,line-style]
      [#:alpha . ,alpha]
      [#:label . ,label]))
  (list (keyword-apply/dict function-interval settings
                            x-axis f
                            x-min* x-max* null)
        (keyword-apply/dict function-interval settings
                            x-axis (invert f)
                            x-min x-max* null)))

(define ((invert f) x)
  (- (f x)))

(define (silverman vals)
  (define iqr (interquartile-range vals))
  (define n (length vals))
  (* 0.9
     (min (stddev vals) (/ iqr 1.34))
     (expt n -0.2)))

(define (interquartile-range vals)
  (- (quantile 3/4 < vals)
     (quantile 1/4 < vals)))

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; example

(parameterize ([plot-y-ticks no-ticks]
               [plot-y-label #f]
               [plot-x-far-ticks no-ticks]
               [plot-x-label "Time (sec)"])
  (plot (list (violin '(0 1 1 2 3 4 4 4 5 6 7 9 10 10 10 11 13))
              (violin '(15 16 17 18 19 20 20 21 23 30)))))