Recursively calculate a list average

Viewed 1577

I have homework to do in OCaml, and one question is about calculating the average of a list. I did that already 1 or 2 years ago in another language, and as I did the first time, I decided to not only sum all elements and divide by the length. The main reason is the fear of a floating-point overflow.

So I found the formula I used the last time on Wikipedia : recursive average formula.

I coded it this way in OCaml :

let average = function
| []    -> raise Empty_list
| hd::l ->
    let rec aux average count = function
        | hd::l -> aux ((average*.(float (count-1))+.hd)/.(float (count))) (count+1) l
        | _     -> average
    in aux hd 1 l
;;

which for me looks like the exact transcription of the formula in OCaml.

But it didn't work, however, after taking a paper, a pen and thinking about it I managed to get it working by replacing the line:

| hd::l -> aux ((average*.(float (count-1))+.hd)/.(float (count))) (count+1) l

with:

| hd::l -> aux ((average*.(float (count))+.hd)/.(float (count+1))) (count+1) l

and it worked.

I told myself the second line is logically the good one to compute the right answer, but I can't understand what was wrong in the first place. Did I translate a biased formula? Or did I miss something when translating?

At this point, it still looks for me like the first line is the transcription of the formula, and the second line is the way to compute the right answer. But I believe there's something I can't understand here. May someone shed light on this for me?

6 Answers

For reference, here is a version of the function that does not overflow with the right time complexity:

let avg l =
  let mu_n' (n,mu_n) x =
    let n' = n + 1 in
    n', mu_n +. (x -. mu_n) /. float n' in
  snd (List.fold_left mu_n' (0,0.) l)

let x = avg [max_float; 1.; 2.; max_float;2.; 3.; max_float; 5.; 6.]
let relative_error = (x -. max_float /. 3.) /. (max_float /. 3.)

val relative_error : float = -1.66533453693773481e-16

But I believe there's something I can't understand here

Nothing is wrong with your logic in general, the formula itself is the source of confusion, I think.

It is quite obvious that (n - 1) multiplier in the dividend must NOT turn into zero during the calculation (otherwise you "discard" the previously accumulated value - which actually happened with your first attempt), and the only way to ensure this is to set n > 0. So, the first equation (the default case) must be indexed by 1, not 0.

So, you have n = 1 for the base case, n = 2 for the next iteration etc. which matches your second (correct) expression and not the first one...

There's a cleaner form of the formula for average which finds the delta between the old average and the new observation, then scales the delta by the sample size to update the average. The base case is that the average of a single observation is that observation. (The average of an empty list is undefined.)

In OCaml:

let rec avg lst =
  match lst with
    | [x]     -> x
    | x::rest -> avg rest +. (x -. avg rest) /. float(List.length lst)
    | []      -> failwith "avg called on empty list!"
;;

The recursive call should only evaluate once since it is pure.

The problem is not in the formula but in the way you use it.

You call aux hd 1 l. So you start with an average of the head of the list and a count of 1. But in the formula you multiply the previous average with count - 1, which is 0 on the first call. So what you do is throw away the head.

Written that way the way to call it is aux 0.0 1 (hd::tl) or aux hd 2 tl.

If you further allow that the average of an empty list is 0.0 you don't even need a pattern matching for the outer function. Going one step further if you make the average and count arguments optional (default to 0.0 and 1 respectively) you don't even need a helper function:

let rec average ?(avg=0.0) ?(count=1) = function
| []     -> avg
| hd::tl -> average
                ~avg:((avg*.(float (count-1))+.hd)/.(float (count)))
                ~count:(count+1)
                tl;;
val average : ?avg:float -> ?count:int -> float list -> float = <fun>

# average [1.;2.;3.];;
- : float = 2.

Why make it so complicated? Why not just calculate the sum and count?

let int_avg lst =
  let rec int_avg_aux cnt sum lst =
    match lst with
    | [] -> (cnt, sum)
    | hd::tl -> int_avg_aux (cnt + 1) (hd + sum) tl in
  int_avg_aux 0 0 lst

let (c, s) = int_avg [1;2;3;4;5;]

let () = Printf.printf "%d %d\n" c s

Now you have the number of elements and the sum of the elements.

I tried your formula in OCaml and I think I got it right:

let avg c lst =
  let rec avg_aux c l =
  match l with
  | [] -> 0.0
  | hd::tl ->
    (((avg_aux (c -. 1.0) tl) *. (c -. 1.0)) +. hd) /. c in
  avg_aux c lst

let lst = [max_float;2.0;max_float;4.0;5.0;6.0]

let ans = avg (float(List.length lst)) lst

let () = Printf.printf "%f\n" ans

Is this what you are looking for?

Related