let chi_squared distrib1 nb_samples1 distrib2 nb_samples2 =
    let length = Array.length distrib1 in
    let degrees_of_freedom = ref length in
    assert (Array.length distrib2 = length);
(*     assert (Array.fold_left (+) 0 distrib1 = nb_samples1); *) (* !! *)
(*     assert (Array.fold_left (+) 0 distrib2 = nb_samples2); (\* time consuming *\) *)
    let normalizer1 = sqrt (float nb_samples2 /. float nb_samples1) in
(*     let normalizer2 = 1. /. normalizer1 in *)
    let sum = ref 0. in
    for i=0 to pred length do
      let x = distrib1.(i) and y = distrib2.(i) in
      if x = 0 && y = 0
(*       if x < 5 && y < 5 *)
      then decr degrees_of_freedom
      else
        let x,y = float x, float y in
(*         let diff = x  -. y in *)
(*         let diff = normalizer1 *. x  -.  normalizer2 *. y in *)
        let diff = normalizer1 *. x  -.  y /. normalizer1 in
        sum := !sum +. diff *. diff /. (x +. y)
    done;
(*     gamma_complement (0.5 *. float !degrees_of_freedom) (0.5 *. !sum) *)
    !sum, !degrees_of_freedom