That's a different algorithm, one with poor numeric properties. For example,
%:((*<:)@#%~(#*+/@:*:)-*:@:(+/)) 1e6+?1e6#0 0.413838 %:((*<:)@#%~(#*+/@:*:)-*:@:(+/)) 1e6+?1e6#0 0.145163 compared with the stddev from stats, which does a direct sum of squared deviance: stddev 1e6+?1e6#0 0.288728 stddev 1e6+?1e6#0 0.28884 Welford's method should have better properties than the sum of squared deviance, but it has a strong sequential-update flavor that I don't know how to do in J. P.S. As should be clear, I dropped a square root from my stddev calculation below. "R.E. Boss" <[email protected]> writes: > ((*<:)@#%~(#*+/@:*:)-*:@:(+/)) > > > R.E. Boss > > >> -----Oorspronkelijk bericht----- >> Van: [email protected] [mailto:programming- >> [email protected]] Namens Johann Hibschman >> Verzonden: dinsdag 13 september 2011 21:23 >> Aan: [email protected] >> Onderwerp: [Jprogramming] Welford's method for standard deviations >> >> I just came across Welford's method for calculating standard deviations, >> and I realized I didn't know how to implement it in J. >> >> See >> >> - http://en.wikipedia.org/wiki/Algorithms_for_calculating_variance >> - http://www.johndcook.com/blog/2008/09/26/comparing-three-methods- >> of-computing-standard-deviation/ >> - http://www.johndcook.com/standard_deviation.html >> >> for the various articles that got me wondering about this. >> >> Here's the algorithm, as stated in the above links. >> >> 1. Initialize m_1 = x_1 and s_1 = 0. >> 2. For subsequent xs, >> m_k = m_{k-1} + (x_k - m_{k-1})/k >> s_k = s_{k-1} + (x_k - m_{k-1})*(x_k - m_k) >> 3. Then the kth running standard deviation is: >> stddev_k = s_k / (k - 1) >> >> I've not tried to switch it from 1-based indexing to 0-based indexing, >> or to use J-like phrasing in the problem, for fear I'd make a stupid >> mistake. >> >> How would I do this in J? >> >> Cheers, >> Johann >> >> ---------------------------------------------------------------------- >> For information about J forums see http://www.jsoftware.com/forums.htm > > ---------------------------------------------------------------------- > For information about J forums see http://www.jsoftware.com/forums.htm ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm
