Re: [fpc-devel] Numeric error in the calculation of Variance

[Werner Pamler]

Am 23.12.2017 um 13:24 schrieb Anton Shepelev:
Hello, all
May  I  ask you if there is any chance of fixing the
error in the Math unit that I reported here:

    https://bugs.freepascal.org/view.php?id=32804

and provided a patch?  Here is a test case:

Program Vartest;
uses Math;

const Size = 1000000;
var dataS: array of Single;
     dataD: array of Double;
     dataE: array of Extended;
     i,n: longint;
begin
   WriteLn('Each run should return a value near unity.');

   WriteLn('Single:');
   SetLength( dataS, Size );
   for n := 1 to 4 do
   begin
     for i := 0 to Size - 1 do
     begin
       dataS[i] := 10000000 + RandG(0,1);
     end;
     WriteLn( Math.Variance( dataS ):5:3 );
   end;

   WriteLn('Double:');
   SetLength( dataD, Size );
   for n := 1 to 4 do
   begin
     for i := 0 to Size - 1 do
     begin
       dataD[i] := 1000000000000000 + RandG(0,1);
     end;
     WriteLn( Math.Variance( dataD ):5:3 );
   end;

   WriteLn('Extended:');
   SetLength( dataE, Size );
   for n := 1 to 4 do
   begin
     for i := 0 to Size - 1 do
     begin
       dataE[i] := 1000000000000000000 + RandG(0,1);
     end;
     WriteLn( Math.Variance( dataE ):5:3 );
   end;
end.

Yes, the formula used is  SS-Sqr(S)/N, and this is the bad one. 
Explained somewhere in the Numerical Recipes, but also in Widepedia 
(https://en.wikipedia.org/wiki/Variance#Formulae_for_the_variance): "it 
is not advised for computer calculations as it suffers from catastrophic 
cancellation <https://en.wikipedia.org/wiki/Catastrophic_cancellation> 
if the two components of the equation are similar in magnitude and 
floating point arithmetic is used."
_______________________________________________
fpc-devel maillist  -  fpc-devel@lists.freepascal.org
http://lists.freepascal.org/cgi-bin/mailman/listinfo/fpc-devel