help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: var(x) and std(x)


From: Paul Kienzle
Subject: Re: var(x) and std(x)
Date: Tue, 12 Nov 2002 13:32:17 -0500
User-agent: Mutt/1.2.5.1i

On Tue, Nov 12, 2002 at 12:17:30PM -0600, Mike Miller wrote:
> On Tue, 12 Nov 2002, Dave Cottingham wrote:
> 
> > I see several people have already explained the n - 1, but I have
> > another concern.  The method cited here is prone to bad roundoff error
> > in cases where sqrt(variance) is much smaller than the mean.  A better
> > method would be
> >
> >    y = sumsq(x - sum(x) / n) / (n - 1);
> 
> 
> This is what it was doing...

std is defined that way (I made sure of that long ago).  var apparently is
not.  Yes, it does make a big difference:

octave:31> q=randn(100,1);
octave:32> x = q + 0;
octave:33> y = (sumsq (x) - sum(x)^2 / n) / (n-1)
y = 0.94445
octave:34> y = sumsq(x - sum(x)/n) / (n-1)
y = 0.94445
octave:35> x=q+1e8;
octave:36> y = (sumsq (x) - sum(x)^2 / n) / (n-1)
y = 11.636
octave:37> y = sumsq(x - sum(x)/n) / (n-1)
y = 0.94445


Index: var.m
===================================================================
RCS file: /cvs/octave/scripts/statistics/base/var.m,v
retrieving revision 1.6
diff -c -r1.6 var.m
*** var.m       2002/09/26 22:48:24     1.6
--- var.m       2002/11/12 18:30:44
***************
*** 40,48 ****
      y = 0;
    elseif ((nr == 1) || (nc == 1))
      n = length (x);
!     y = (sumsq (x) - abs(sum(x))^2 / n) / (n - 1);
    else
!     y = (sumsq (x) - abs(sum(x)).^2 / nr) / (nr - 1);
    endif
  
  endfunction
--- 40,48 ----
      y = 0;
    elseif ((nr == 1) || (nc == 1))
      n = length (x);
!     y = sumsq (x - sum (x) / n) / (n - 1);
    else
!     y = sumsq (x - ones (nr, 1) * (sum (x) / nr) ) / (nr - 1);
    endif
  
  endfunction



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------



reply via email to

[Prev in Thread] Current Thread [Next in Thread]