;; For background, see http://www.johndcook.com/standard_deviation.html

(define (var-welford vec)
  "Calculate sample variance of vector by Welford's method."
  (let loop ((n (f64vector-length vec))
	     (k 1) 
	     (m (f64vector-ref vec 0))
	     (s 0))
    (if (= k n)
	(/ s (sub1 n))
	(let* ((xk (f64vector-ref vec k))
	       (mk (+ m (/ (- xk m) (add1 k))))
	       (sk (+ s (* (- xk m) (- xk mk)))))
	  (loop n (add1 k) mk sk)))))