wisp
 
(Arne Babenhauserheide)
2017-03-21: Plot fit results for datastructure scaling behavior

Plot fit results for datastructure scaling behavior

diff --git a/examples/benchmark.w b/examples/benchmark.w
--- a/examples/benchmark.w
+++ b/examples/benchmark.w
@@ -55,20 +55,12 @@ define-syntax benchmark
         #' benchmark-fun
          . (lambda () thunk) args ...
 
-;; TODO: Use fit to different mappings.
-define : mismatch-to-const-N-m timing-list
-  define : N-m x
-    define : const y
-             car : cdr x
-    map const : car x
-  map N-m timing-list
-
-define : mismatch-to-linear-N-m timing-list
-  define : N-m x 
-    define : linear y
-       / (car (cdr x)) y
-     map linear : car x
-  map N-m timing-list
+define : logiota steps start stepsize
+    . "Create numbers evenly spread in log space"
+    let*
+        : logstart : log (+ start 1)
+          logstep : / (- (log (+ start (* stepsize (- steps 1)))) logstart) (- steps 1)
+        map inexact->exact : map round : map exp : iota steps logstart logstep 
 
 define : benchmark-list-append
   . "Test (append a b) with lists of different lengths."
@@ -79,21 +71,21 @@ define : benchmark-list-append
            let : (N (list-ref x 0)) (m (list-ref x 1))
                benchmark (append a b) :let ((a (iota N))(b (iota m)))
          . param-list
-  let : (steps 30)
+  let : (steps 50)
     concatenate
       list 
-        let : (param-list (zip (iota steps 1 1000) (iota steps 1 0)))
-               bench-append param-list
-        let : (param-list (zip (iota steps 1 0) (iota steps 1 1000)))
-               bench-append param-list
-        let : (param-list (zip (iota steps 1 1000) (iota steps 1 0)))
+        let : (param-list (zip (logiota steps 1 100) (logiota steps 1 0)))
                bench-append param-list
-        let : (param-list (zip (iota steps 1 0) (iota steps 1 1000)))
-               bench-append param-list
-        let : (param-list (zip (iota steps 1 1000) (iota steps 100000 0)))
-               bench-append param-list
-        let : (param-list (zip (iota steps 100000 0) (iota steps 1 1000)))
-               bench-append param-list
+        ;; let : (param-list (zip (logiota steps 100 0) (logiota steps 1 100)))
+        ;;        bench-append param-list
+        ;; let : (param-list (zip (logiota steps 1 1000) (logiota steps 1 0)))
+        ;;        bench-append param-list
+        ;; let : (param-list (zip (logiota steps 1 0) (logiota steps 1 1000)))
+        ;;        bench-append param-list
+        ;; let : (param-list (zip (logiota steps 1 1000) (logiota steps 100000 0)))
+        ;;        bench-append param-list
+        ;; let : (param-list (zip (logiota steps 100000 0) (logiota steps 1 1000)))
+        ;;        bench-append param-list
 
 ;; stddev from rosetta code: http://rosettacode.org/wiki/Standard_deviation#Scheme
 define : stddev nums
@@ -128,24 +120,24 @@ x are parameters to be optimized, pos is
            +
              list-ref x 0 ; constant value
              ;; pure N
-             * (list-ref x 1) : log : + 1 N ; avoid breakage at pos 0
-             ; * (list-ref x 2) : sqrt N
+             * (list-ref x 1) : log (+ 1 N) ; avoid breakage at pos 0
+             * (list-ref x 2) : sqrt N
              * (list-ref x 3) N
-             ; * (list-ref x 4) : expt N 2
-             ; * (list-ref x 5) : expt N 3
+             * (list-ref x 4) : * N : log (+ 1 N)
+             * (list-ref x 5) : expt N 2
              ;; pure m
-             * (list-ref x 6) : log : + 1 m ; avoid breakage at pos 0
-             ; * (list-ref x 7) : sqrt m
+             * (list-ref x 6) : log (+ 1 m) ; avoid breakage at pos 0
+             * (list-ref x 7) : sqrt m
              * (list-ref x 8) m
-             ; * (list-ref x 9) : expt m 2
-             ; * (list-ref x 10) : expt m 3
+             * (list-ref x 9) : * m : log (+ 1 m)
+             * (list-ref x 10) : expt m 2
              ;; mixed terms
-             * (list-ref x 11) : log : + 1 N m
-             * (list-ref x 12) : * N (log (+ 1 m))
-             * (list-ref x 13) : * m (log (+ 1 N))
-             ; * (list-ref x 14) : * N m
-             ; * (list-ref x 15) : * (expt N 2) m
-             ; * (list-ref x 16) : * (expt m 2) N
+             * (list-ref x 11) : log (+ 1 N m)
+             * (list-ref x 12) : * N : log (+ 1 m)
+             * (list-ref x 13) : * m : log (+ 1 N)
+             * (list-ref x 14) : * N m
+             * (list-ref x 15) : * (expt N 2) m
+             * (list-ref x 16) : * (expt m 2) N
 
 
 define : interleave lx lz
@@ -157,13 +149,21 @@ define : interleave lx lz
 
 
 define : print-fit x σ
-    let : : msg "
-~,1,,,,,'ee±~,1,,,,,'ee + ~,1,,,,,'ee±~,1,,,,,'ee log(N) + ~,1,,,,,'ee±~,1,,,,,'ee sqrt(N) + ~,1,,,,,'ee±~,1,,,,,'ee N + ~,1,,,,,'ee±~,1,,,,,'ee N^2 + ~,1,,,,,'ee±~,1,,,,,'ee N^3 +
-~,1,,,,,'ee±~,1,,,,,'ee log(m) + ~,1,,,,,'ee±~,1,,,,,'ee sqrt(m) + ~,1,,,,,'ee±~,1,,,,,'ee m + ~,1,,,,,'ee±~,1,,,,,'ee m^2 + ~,1,,,,,'ee±~,1,,,,,'ee m^3 +
-~,1,,,,,'ee±~,1,,,,,'ee log(N + m) + ~,1,,,,,'ee±~,1,,,,,'ee N log(m) + ~,1,,,,,'ee±~,1,,,,,'ee m log(N)+ ~,1,,,,,'ee±~,1,,,,,'ee N m + ~,1,,,,,'ee±~,1,,,,,'ee N^2 m + ~,1,,,,,'ee±~,1,,,,,'ee m^2 N
-"
-      apply format 
-        append (list #t msg) (interleave x σ)
+    . "Print the big-O parameters which are larger than σ (their standard deviation)."
+    let : : number-format "~,1,,,,,'ee±~,1,,,,,'ee"
+      let big-O
+        : names : list "" "log(N)" "sqrt(N)" "N log(N)" "N^2" "log(m)" "sqrt(m)" "m" "m log(m)" "m^2" "log(N + m)" "N log(m)" "m log(N)" "N m" "N^2 m" "m^2 N"
+          x x
+          σ σ
+        cond
+          : or (null? names) (null? x) (null? σ)
+            newline
+          : > (abs (car x)) (car σ)
+            format #t : string-append number-format " " (car names) "  "
+                      . (car x) (car σ)
+            big-O (cdr names) (cdr x) (cdr σ)
+          else
+            big-O (cdr names) (cdr x) (cdr σ)
 
 
 define : flatten li
@@ -180,15 +180,19 @@ define : main args
    ;;     list mismatch-to-const-N-m mismatch-to-linear-N-m
    let*
       : bench : benchmark-list-append ;; benchmark results
-        ; fitting to cost estimates
-        ensemble-member-count 32
+        ;; fitting to cost estimates
+        ensemble-member-count 128
         ensemble-member-plot-skip 4
-        x^b : list-ec (: i 17) (car (cdr (car bench))) ; inital guess: constant starting at the first result
-        x^b-std : list-ec (: i 17) (car (cdr (car bench))) ; inital guess: 100% uncertainty
+        y_0 : apply min : map car : map cdr bench
+        y_m : apply max : map car : map cdr bench
+        nb : apply max : interleave (map car (map car bench)) (map car (map cdr (map car bench)))
+        ;; "const" "log(N)" "sqrt(N)" "N" "N^2" "N^3" "log(m)" "sqrt(m)" "m" "m^2" "m^3" "log(N + m)" "N log(m)" "m log(N)" "N m" "N^2 m" "m^2 N"
+        x^b : list y_0 (/ y_m (log nb)) (/ y_m (sqrt nb)) (/ y_m nb) (/ y_m nb nb) (/ y_m nb nb nb) (/ y_m (log nb)) (/ y_m (sqrt nb)) (/ y_m nb) (/ y_m nb nb) (/ y_m nb nb nb) (/ y_m nb nb) (/ y_m nb nb) (/ y_m nb nb nb) (/ y_m nb nb nb) (/ y_m nb nb nb nb) (/ y_m nb nb nb nb) ; inital guess: constant starting at the first result
+        x^b-std : list-ec (: i x^b) i ; inital guess: 100% uncertainty
         P : make-covariance-matrix-with-offdiagonals-using-stds x^b-std
         y⁰-pos : map car bench
         y⁰ : append-map cdr bench
-        y⁰-std : stddev y⁰
+        y⁰-std : list-ref (sort y⁰ <) : round : / (length y⁰) 16 ; lower octile median
         R : make-covariance-matrix-with-offdiagonals-using-stds : list-ec (: i (length bench)) y⁰-std
         optimized : EnSRF H x^b P y⁰ R y⁰-pos ensemble-member-count
         x-opt : list-ref optimized 0
@@ -224,28 +228,33 @@ define : main args
         format port "yinitstds = [float(i) for i in '~A'[1:-1].split(' ')]\n" y^b-stds
         format port "yopt = [float(i) for i in '~A'[1:-1].split(' ')]\n" : list-ec (: i y⁰-pos) : H x-opt i
         format port "yoptstds = [float(i) for i in '~A'[1:-1].split(' ')]\n" y-stds
-        format port "pl.errorbar(*zip(*sorted(zip(ypos1, yinit))), yerr=zip(*sorted(zip(ypos1, yinitstds)))[1], label='prior vs N')\n"
-        format port "pl.errorbar(*zip(*sorted(zip(ypos1, yopt))), yerr=zip(*sorted(zip(ypos1, yoptstds)))[1], label='optimized vs N')\n"
+        ;; format port "pl.errorbar(*zip(*sorted(zip(ypos1, yinit))), yerr=zip(*sorted(zip(ypos1, yinitstds)))[1], label='prior vs N')\n"
+        format port "pl.errorbar(*zip(*sorted(zip(ypos1, yopt))), yerr=zip(*sorted(zip(ypos1, yoptstds)))[1], marker='+', mew=2, ms=10, linewidth=0.1, label='optimized vs N')\n"
         format port "eb=pl.errorbar(*zip(*sorted(zip(ypos1, y0))), yerr=yerr, alpha=0.6, marker='x', mew=2, ms=10, linewidth=0, label='measurements vs N')\neb[-1][0].set_linewidth(1)\n"
-        format port "pl.errorbar(*zip(*sorted(zip(ypos2, yinit))), yerr=zip(*sorted(zip(ypos2, yinitstds)))[1], label='prior vs. m')\n"
-        format port "pl.errorbar(*zip(*sorted(zip(ypos2, yopt))), yerr=zip(*sorted(zip(ypos2, yoptstds)))[1], label='optimized vs. m')\n"
+        ;; format port "pl.errorbar(*zip(*sorted(zip(ypos2, yinit))), yerr=zip(*sorted(zip(ypos2, yinitstds)))[1], label='prior vs. m')\n"
+        format port "pl.errorbar(*zip(*sorted(zip(ypos2, yopt))), yerr=zip(*sorted(zip(ypos2, yoptstds)))[1], marker='+', mew=2, ms=10, linewidth=0.1, label='optimized vs. m')\n"
         format port "eb=pl.errorbar(*zip(*sorted(zip(ypos2, y0))), yerr=yerr, alpha=0.6, marker='x', mew=2, ms=10, linewidth=0, label='measurements vs. m')\neb[-1][0].set_linewidth(1)\n"
-        list-ec (: step 0 (length x^steps) 16)
+        format port "pl.plot(sorted(ypos1), pl.log(sorted(ypos1))*(max(y0) / pl.log(max(ypos1))), label='log(N)')\n"
+        format port "pl.plot(sorted(ypos1), pl.sqrt(sorted(ypos1))*(max(y0) / pl.sqrt(max(ypos1))), label='sqrt(N)')\n"
+        format port "pl.plot(sorted(ypos1), pl.multiply(sorted(ypos1), max(y0) / max(ypos1)), label='N')\n"
+        list-ec (: step 0 (length x^steps) 4)
              let : : members : list-ref x^steps (- (length x^steps) step 1)
                 list-ec (: member-idx 0 (length members) ensemble-member-plot-skip) ; reversed
                    let : : member : list-ref members member-idx
                      format port "paired = pl.get_cmap('Paired')
 cNorm = mpl.colors.Normalize(vmin=~A, vmax=~A)
 scalarMap = mpl.cm.ScalarMappable(norm=cNorm, cmap=paired)\n" 0 (length member)
-                     list-ec (: param-idx 0 (length member) 16) ; step = 16
-                        ; plot parameter 0
-                        format port "pl.plot(~A, ~A, marker='.', color=scalarMap.to_rgba(~A), linewidth=0, label='', alpha=0.6, zorder=-1)\n" (/ step 1) (+ 80 (* (/ (apply + y-opt) (length y-opt)) (list-ref member param-idx))) param-idx
+                     list-ec (: param-idx 0 (length member) 4) ; step = 4
+                        ;; plot parameter 0
+                        let : (offset (/ (apply max (append y⁰ y-opt)) 2)) (spreading (/ (apply max (append y⁰ y-opt)) (- (apply max member) (apply min member))))
+                            format port "pl.plot(~A, ~A, marker='.', color=scalarMap.to_rgba(~A), linewidth=0, label='', alpha=0.6, zorder=-1)\n" 
+                                        . (/ step 1) (+ offset (* spreading (list-ref member param-idx))) param-idx
         format port "pl.legend(loc='upper right')\n"
         format port "pl.xlabel('position [arbitrary units]')\n"
         format port "pl.ylabel('value [arbitrary units]')\n"
-        format port "pl.title('ensemble optimization results')\n"
+        format port "pl.title('~A')\n" "Operation scaling behaviour"
         format port "pl.xscale('log')\n"
-        format port "pl.yscale('log')\n"
+        ;; format port "pl.yscale('log')\n"
         format port "pl.show()\n"
         format port "exit()\n"
         close-pipe port