Thanks to everyone who responded. Mike, I like your technique of returning output in the form of input for a future iteration.
On Thu, Sep 20, 2018 at 6:22 PM 'Mike Day' via Programming < [email protected]> wrote: > It's likely my old stuff was in APL. > > But for now, I've worked up the following code which should work for > several > groups of data, stored in three columns containing numbers, means and > unbiased > standard deviations. > > (Later thought while drafting - I think it should work for biased estimates > too, if we replace n - 1, below, by n - UNBIASED where global UNBIASED is > 1 or 0 respectively. Recode accordingly if you wish!) > > NB. retrieve sum of sqs of x given n, m, s for one set > sx2 =: 3 : 0 > 'n m s' =. y > (n - 1 0) (+/ . *) *: s, m > ) > > NB. sd given n, m, sx2 > sd =: 3 : 0 > 'n m s2' =. y NB. number, mean, sum of squared values > %: (n - 1) %~ s2 - n (* *:) m > ) > > NB. Given n1, m1, s1, n2, m2, s2, etc, in 3 columns, > NB. what is overall n,m,s for two sets combined? > nmscomb =: 3 : 0 > 'n m s' =. |:y > sumxx =. +/sx2"1 y NB. overall sum of squares of values > nmall =. ({:, %/) n +/ . * m,. 1 NB. overall number & mean > (}:, sd) nmall, sumxx NB. append overall standard > deviation > : > nmscomb x,:y NB. for two sets in lh & rh args > ) > > 5 51.4 23.6072 nmscomb 12 70.4167 26.3972 NB. Devon's "data" as 2 sets > 17 64.8236 26.4226 > > 2 3 $ 5 51.4 23.6072 12 70.4167 26.3972 NB. as table of with rows > 5 51.4 23.6072 > 12 70.4167 26.3972 > > nmscomb 2 3 $ 5 51.4 23.6072 12 70.4167 26.3972 NB. process table > 17 64.8236 26.4226 > > Any use? > > Mike > > > On 20/09/2018 18:01, 'Mike Day' via Programming wrote: > > I was going to have a look back at my old workings from when I was > > actually > > involved professionally in thinking about meta-stats, but you've > > evidently > > found something useful while we were having supper! > > > > Let me know if you think it might still be useful, always assuming I > > can find > > anything relevant! > > > > Cheers, > > > > Mike > > > > [NB no snipping of back-thread!] > > > > On 20/09/2018 17:20, Devon McCormick wrote: > >> OK - it turns out that there's a very simple solution, courtesy of this > >> site: > >> > https://www.mathworks.com/matlabcentral/fileexchange/35605-overall-mean-standard-deviation-of-groups-of-observations?w.mathworks.com > >> > >> . > >> > >> combineSDsScalars=: 3 : 0 > >> 'grpn grpmean grpstd'=. y > >> alln=. +/grpn NB. Total # observations > >> allmean=. grpmean +/ . * grpn%alln NB. Weighted mean of all > >> ESS=. (*:grpstd) +/ . * <:grpn NB. Total error sum of squares > >> GSS=. grpn +/ . * *:grpmean-allmean NB. Total group sum of squares > >> %:(ESS+GSS)%<:alln NB. Total standard deviation > >> ) > >> > >> nn=. 23 81 52 67 34; 97 72 53 22 85 39 97 94 82 35 79 90 > >> <"1|:(#,mean,stddev)&>nn > >> +----+------------+---------------+ > >> |5 12|51.4 70.4167|23.6072 26.3972| > >> +----+------------+---------------+ > >> combineSDsScalars <"1|:(#,mean,stddev)&>nn > >> 26.4226 > >> stddev ;nn NB. Correct answer > >> 26.4226 > >> nn2=. 59 52 90 62 0; 64 29 14 30 92 98 82 12 54 30; 84 55 35 7 8 > >> 23 66 > >> 60 99 62 93 91 35 41 27 > >> stddev ;nn2 NB. Correct answer > >> 30.2466 > >> combineSDsScalars <"1|:(#,mean,stddev)&>nn2 > >> 30.2466 > >> > >> > >> > >> > >> > >> On Thu, Sep 20, 2018 at 11:40 AM Raul Miller <[email protected]> > >> wrote: > >> > >>> I think the approach should extend to multiple series if you build a > >>> combine routine which works in a reduce context. > >>> > >>> In other words, you need to build all the summary statistics when > >>> combining two series. Or, something like this (warning: untested > >>> code): > >>> > >>> combineSum=: +&{. > >>> > >>> combineMean=: +/&(*/)&}: % combineSum > >>> > >>> combine=: combineSum, combineMean, combineSDsS@,@,. > >>> > >>> I *think* combineSDsS would be like combineSDsScalars, but without the > >>> m argument (which doesn't seem to be relevant). If that's right, you'd > >>> use it like this: > >>> > >>> combine&;/summary1;summary2;summary3 > >>> > >>> I hope this helps, > >>> > >>> -- > >>> Raul > >>> > >>> On Thu, Sep 20, 2018 at 11:18 AM Devon McCormick <[email protected]> > >>> wrote: > >>>> I work with a system that handles a lot of complicated problems very > >>> simply > >>>> but I've recently run into a difficult problem, the solution of which > >>> would > >>>> be how do we re-construct the overall standard deviation of a > >>>> number of > >>>> series, to which we would rather not refer to in detail, when we have > >>> only > >>>> the sample size, the mean, and the standard deviation for each group. > >>>> > >>>> How do we combine these? > >>>> > >>>> There's a lot of stuff on the web but most of what I've found is > >>>> wrong - > >>>> unless I'm wrong and the numerous different answers are all right - > >>>> but > >>> for > >>>> this reference - almost: > >>>> > >>> > http://atozmath.com/CONM/Ch2_CombinedSD.aspx?q1=1.1%605%2c12%6051.4%2c70.4167%6023.6072%2c26.3972%60SD2#PrevPart > >>> > >>> . > >>>> I say "almost" because, on the website, the formula is correct but the > >>>> arithmetic is wrong! > >>>> > >>>> The site shows an example worked out in painstaking detail - given the > >>>> summary stats for this example - but the final answer is shown as > >>>> 26.1952 > >>>> when the correct answer is 26.4226. > >>>> > >>>> However, if I properly reproduce the formula and do my own math, I get > >>> the > >>>> correct answer, to whit - here's the formula, as given, with the math, > >>> and > >>>> following is my own work in J, reproducing the formulas but getting > >>>> the > >>>> math right. > >>>> Combined Mean = 64.8236 > >>>> Combined Standard deviation : > >>>> σ12=√(N1-1)⋅σ21+(N2-1)⋅σ22+N1⋅N2N1+N2⋅(ˉx21+ˉx22-2ˉx1ˉx2)N1+N2-1 > >>>> ]nn=. (5?@$100);12?@$100 > >>>> +--------------+-----------------------------------+ > >>>> |23 81 52 67 34|97 72 53 22 85 39 97 94 82 35 79 90| > >>>> +--------------+-----------------------------------+ > >>>> stddev ;nn NB. Correct answer > >>>> 26.4226 > >>>> (*:@:stddev)&>nn > >>>> 557.3 696.811 > >>>> var&>nn > >>>> 557.3 696.811 > >>>> > >>>> (#,mean,stddev)&>nn > >>>> 5 51.4 23.6072 > >>>> 12 70.4167 26.3972 > >>>> > >>>> NB. Following based on formula here: > >>>> > >>> > http://atozmath.com/CONM/Ch2_CombinedSD.aspx?q1=1.1%605%2c12%6051.4%2c70.4167%6023.6072%2c26.3972%60SD2#PrevPart > >>> > >>>> combineSDs=: [: %: ((([: <: #&>) +/ .* var&>) + ([: (*/ % +/) > >>>> #&>) * > >>> ([: > >>>> +/ [: *: mean&>) - [: +: [: */ mean&>) % [: <: [: +/ #&> > >>>> combineSDs nn NB. This works based on the mean, > >>>> variance, and > >>>> counts of nn. > >>>> 26.4226 > >>>> > >>>> NB. This version makes it clear we can do this using only the > >>>> summary > >>>> statistics: > >>>> 'n1 n2'=. #&>nn [ 's21 s22'=. var&>nn [ 'm1 m2'=. mean&>nn [ > >>>> m=. mean > >>> ;nn > >>>> > >>> > %:(((n1-1)*s21)+((n2-1)*s22)+((n1*n2)%n1+n2)*(*:m1)+(*:m2)-+:m1*m2)%<:n1+n2 > >>> > >>>> 26.4226 > >>>> combineSDsScalars=: 3 : 0 > >>>> 'n1 n2 m1 m2 m s21 s22'=. y > >>>> > >>>> > >>> > %:(((n1-1)*s21)+((n2-1)*s22)+((n1*n2)%n1+n2)*(*:m1)+(*:m2)-+:m1*m2)%<:n1+n2 > >>> > >>>> ) > >>>> combineSDsScalars n1,n2,m1,m2,m,s21,s22 > >>>> 26.4226 > >>>> > >>>> But how do I extend this for more than two series? > >>>> > >>>> -- > >>>> > >>>> Devon McCormick, CFA > >>>> > >>>> Quantitative Consultant > >>>> ---------------------------------------------------------------------- > >>>> For information about J forums see > http://www.jsoftware.com/forums.htm > >>> ---------------------------------------------------------------------- > >>> For information about J forums see http://www.jsoftware.com/forums.htm > >> > >> > > > > > > --- > > This email has been checked for viruses by Avast antivirus software. > > https://www.avast.com/antivirus > > > > ---------------------------------------------------------------------- > > For information about J forums see http://www.jsoftware.com/forums.htm > > ---------------------------------------------------------------------- > For information about J forums see http://www.jsoftware.com/forums.htm -- Devon McCormick, CFA Quantitative Consultant ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm
