Weighted standard deviation calculation

[Jesse Tielens]

Dear statalisters,

I hope someone here is able to help me find out how to calculate a weighted standard deviation. I have included a datasample by Dataex at the end of this post, as it's rather long.

I've got panel data on stock performance, firms are identified with an id and have observations for each year. I've calculated Abnormal Returns (AR) and sorted stocks from low to high in quintiles by the variable CI (CI_quintile). My aim is to show the abnormal returns to an investment strategy that invests into all the stocks listed in a CI_quintile. The problem is, if I invest a single Dollar into the portfolio this dollar is divded over all these stocks according to their portfolio weight (varname = weight), because a large company makes up a larger share of the portfolio of course.

An individual stock's weighted Abnormal Return is weight * AR = wAR.

When I want to know the population mean of abnormal returns by group, I can do:

 tabstat AR, by(CI_quintile) statistics(count mean semean median max min) columns(statistics)

This will show me the unweighted mean, median and standard deviation for each CI_quintile. If I want to look at the outcome of a portfolio strategy, I use the following code to calculate the yearly weighted return of an investment strategy into each CI_quintile and create the table below:

  
 gen wAR = weight * AR
bysort year CI_quintile: egen weighted_total_AR = sum(wAR)
tabstat weighted_total_AR, by(CI_quintile) statistics (mean semean median max min)  
//Because weighted_total_AR is an average, but egen generates observations for all. I just want to look at the first occurence of weighted_total in each year/group.
bysort year CI_quintile: egen first = _n tabstat weighted_total_AR if first ==1, by(CI_quintile) statistics (mean semean median max min)

But, the standard deviation displayed in the second column is not the right one. This is calculated using the variance between the mean returns of the portfolios. What I'm actually interested in is the weighted standard deviation of each portfolio. So that Xi - X is the difference between the mean abnormal return in a portfolio and that observation's abnormal return. This way, I can do some predictions about whether the observed abnormal returns are statistically significant, etc. I can calculate this for each year's portfolio by hand, doing this:

So in words: For each category I'd need to multiply the variance with its weight, and divide the total variance by the sum of all weights. Because in my case, M is equal to N.
Doing it by hand is not an option for me, as I will need to do it for a lot of portfolios.

Is the only option in Stata to calculate this weighted standard deviation by typing the formulas in by hand?
I'm hoping someone can help me find a quick way or trick to do this! I've searched the help pages for additional options in calculating standard deviations, etc. but all I found was that the 'Svy: mean' command has additional options, but those are designed for survey research.

Dataexample:

* Example generated by -dataex-. To install: ssc install dataex
clear
input float id int year float(CI_quintile weight AR wAR)
7458 1988 1 .003308191 .04755842 1.5733234e-06
7348 1988 1 .01979438 27.21212 .00538647
6821 1988 1 .000877853 19.46686 .0001708904
7031 1988 1 .004767775 -64.56488 -.003078308
8037 1988 1 .0007281441 -20.45289 -.00014892654
3976 1988 1 .002622888 19.17756 .0005030059
6948 1988 1 .06607924 32.777126 .021658873
8446 1988 1 .0042028353 -63.46425 -.002667298
5475 1988 1 .004709554 -14.628477 -.0006889361
5418 1988 1 .010878066 .01212345 1.318797e-06
6138 1988 1 .0040314705 27.215406 .0010971811
6555 1988 1 .00775756 11.361523 .0008813769
5717 1988 1 .018076962 -16.627876 -.003005815
5496 1988 1 .03864232 -6.418726 -.002480345
7254 1988 1 .003252011 484.86795 .01576796
7145 1988 1 .07200812 -3.9270766 -.002827814
7294 1988 1 .003709455 -40.03103 -.001484933
6462 1988 1 .018823624 -1.3978765 -.00026313102
5870 1988 1 .008651262 -13.30874 -.001151374
6340 1988 1 .008757503 4.581522 .0004012269
7923 1988 1 .0027628676 -59.12244 -.001633475
7154 1988 1 .001613059 -40.30289 -.0006501094
6665 1988 1 .009262182 -24.833794 -.0023001512
5974 1988 1 .09308114 .786275 .0007318737
6946 1988 1 .032222416 47.52064 .015312297
5604 1988 1 .0037643795 -3.112319 -.0001171595
8504 1988 1 .004078392 15.433672 .0006294457
6836 1988 1 .0020328404 45.48768 .000924692
7217 1988 1 .013693974 32.054104 .004389481
6607 1988 1 .002903632 20.31768 .0005899506
6495 1988 1 .003260328 8.757681 .00028552912
6958 1988 1 .013698054 18.946796 .0025953425
5523 1988 1 .04979282 3.9112735 .0019475333
7184 1988 1 .003675088 65.49897 .002407145
5446 1988 1 .0123105 -19.09788 -.002351044
6859 1988 1 .003496818 -35.66232 -.0012470464
5397 1988 1 .0020494745 -46.04244 -.0009436281
5382 1988 1 .0031159546 45.20768 .001408651
5781 1988 1 .0020171476 -68.512436 -.001381997
7149 1988 1 .013541284 -3.033204 -.0004107348
6056 1988 1 .002672791 -43.60244 -.001165402
5655 1988 1 .05581209 -18.313065 -.010220904
5945 1988 2 .00463354 -1.204877 -.00005582845
6923 1988 2 .04139749 14.199373 .005878184
7144 1988 2 .0020478633 -57.57103 -.001178976
6309 1988 2 .0012699616 2.237681 .00002841769
7311 1988 2 .022528274 33.696278 .00759119
7204 1988 2 .0028297885 -48.55125 -.0013738978
5491 1988 2 .05930467 -13.633725 -.008085435
5824 1988 2 .002104283 -44.38244 -.0009339322
6362 1988 2 .0011264324 -8.682319 -.00009780045
6944 1988 2 .04139796 -7.973727 -.00330096
3869 1988 2 .0022640927 174.80896 .003957837
6647 1988 2 .03457022 46.43628 .016053123
6817 1988 2 .009060114 66.086205 .005987485
6921 1988 2 .02863818 70.34693 .02014608
5417 1988 2 .016566856 -14.123064 -.002339748
8473 1988 2 .00104213 -18.272894 -.0001904273
5412 1988 2 .004920411 -5.548478 -.0002730079
7151 1988 2 .0041635614 22.27852 .00092758
6971 1988 2 .001034832 -25.23244 -.00026111337
6205 1988 2 .0029855745 -24.16848 -.000721568
5461 1988 2 .00498731 -10.508478 -.00052409037
1571 1988 2 .00548985 -43.0032 -.002360811
8454 1988 2 .00383496 8.991261 .0003448113
7312 1988 2 .004747315 -10.678477 -.000506941
5399 1988 2 .02433521 65.471275 .015932571
6826 1988 2 .00601999 -31.88619 -.0019195458
7416 1988 2 .002490614 134.89755 .003359778
5613 1988 2 .0021592996 2.457681 .00005306869
6789 1988 2 .001412181 180.52898 .002549396
6012 1988 2 .006848043 -13.387877 -.0009168075
6913 1988 3 .04084988 31.08292 .012697336
6793 1988 3 .0042613777 54.11652 .0023061095
5769 1988 3 .0044555734 -26.13619 -.0011645172
8524 1988 3 .0017429963 -10.454877 -.0001822281
6201 1988 3 .05990786 -2.725325 -.0016326838
8097 1988 3 .004850253 -25.023205 -.001213689
5485 1988 3 .0787017 -29.518726 -.02323174
7202 1988 3 .003989745 -36.293793 -.0014480298
6800 1988 3 .0079699615 17.086206 .001361764
5756 1988 3 .00467556 -24.70619 -.0011551529
6949 1988 3 .008540733 9.826796 .0008392804
6035 1988 3 .023817774 -26.364063 -.006279333
7427 1988 3 .0014571343 -30.97289 -.0004513166
7628 1988 3 .0019216284 129.47617 .002488051
7166 1988 3 .0007566451 2.0971072 .00001586766
7730 1988 3 .001180102 9.018971 .00010643308
6017 1988 3 .06806776 2.056275 .0013996603
5827 1988 3 .0040634973 -1.8178765 -.00007386936
7428 1988 3 .0021379942 -68.524254 -.0014650446
6939 1988 3 .00689385 -.673204 -.00004640967
7362 1988 4 .005687082 10.156796 .00057762535
7756 1988 4 .0012734267 37.337105 .0004754607
5708 1988 4 .013533988 -23.97708 -.003245055
5883 1988 4 .005295252 -50.66379 -.0026827755
7554 1988 4 .005735876 52.59091 .003016549
8444 1988 4 .0030944906 4.485123 .00013879171
7642 1988 4 .004470479 18.65236 .0008338499
7147 1988 4 .0007827005 -5.032893 -.00003939248
7289 1988 4 .011478996 29.783806 .003418882
7251 1988 4 .002458345 68.115746 .00167452
5654 1988 4 .003510936 -30.404877 -.0010674957
7457 1988 4 .0019048592 54.96897 .0010470815
6968 1988 4 .011513387 -36.653202 -.004220025
8363 1988 4 .004779474 22.485746 .0010747005
6037 1988 4 .015742868 28.75292 .0045265346
6558 1988 4 .0015354616 -46.53232 -.0007144859
6918 1988 4 .03322527 19.520636 .006485784
7484 1988 4 .0014729643 -16.942442 -.00024955612
5632 1988 4 .012606042 -8.817077 -.0011114845
7142 1988 4 .0021627042 18.535107 .0004008595
5778 1988 4 .0014651957 -22.92244 -.0003358586
5483 1988 4 .021556865 -24.233727 -.005224032
6000 1988 4 .0014881522 -17.212318 -.00025614552
5885 1988 4 .016856385 -4.5930634 -.0007742245
4063 1988 4 .0015768355 -22.65103 -.00035716945
7650 1988 4 .004959285 81.96282 .0040647704
5400 1988 4 .005623101 3.056206 .00017185355
7162 1988 4 .0018952576 -6.292893 -.00011926654
6950 1988 4 .018039556 14.067924 .002537791
7258 1988 4 .003267144 5.416422 .0001769623
6908 1988 4 .0034998504 2.175123 .00007612605
7238 1988 5 .006276189 -50.9338 -.003196701
7888 1988 5 .002253436 -4.7310286 -.0001066107
5649 1988 5 .06367889 -22.983725 -.014635782
5616 1988 5 .01350955 10.573808 .0014284737
6176 1988 5 .04494316 23.071924 .010369252
6815 1988 5 .0026754425 11.98006 .0003205196
5383 1988 5 .05657749 4.3412733 .002456183
6206 1988 5 .007973183 -8.053204 -.0006420967
5518 1988 5 .0016337364 -41.03232 -.00067036
5626 1988 5 .02612137 -42.73373 -.011162634
6092 1988 5 .002454278 -20.61244 -.0005058866
5748 1988 5 .017545564 34.69292 .006087069
5944 1988 5 .012790716 2.6538076 .000339441
7266 1988 5 .013038692 14.093807 .001837648
5405 1988 5 .002873613 67.12152 .0019288127
7277 1988 5 .0226874 131.61694 .029860463
7672 1988 5 .0014592282 22.51768 .0003285844
7227 1988 5 .005988808 -20.276304 -.001214309
8509 1988 5 .003483581 137.14513 .004777561
6674 1988 5 .02666846 -25.893724 -.006905457
6810 1988 5 .008821299 33.444523 .0029502416
7169 1989 1 .026088715 44.665 .011652525
6304 1989 1 .00771959 -9.9 -.0007642394
6000 1989 1 .0017870337 -8.860001 -.0001583312
5461 1989 1 .009031327 -27.83375 -.002513757
7458 1989 1 .0038669384 20.1116 .0007777033
6816 1989 1 .010155418 -24.95973 -.002534765
7181 1989 1 .0023925565 -17.099154 -.0004091069
5652 1989 1 .01419654 -14.320545 -.002033022
7416 1989 1 .004921658 .3537169 .000017408738
6872 1989 1 .019402415 6.389994 .001239813
5432 1989 1 .0026918836 17.859463 .000480756
5546 1989 1 .0038483934 -18.571756 -.0007147143
5894 1989 1 .0014853714 -22.564 -.0003351592
6926 1989 1 .005381981 -27.289904 -.0014687375
4063 1989 1 .003162235 2.6484454 .00008375008
7154 1989 1 .0014550128 2.8294635 .00004116906
7293 1989 1 .007206242 28.92644 .0020845092
6651 1989 1 .004246764 -19.09472 -.0008109077
3706 1989 1 .003678743 -24.972 -.0009186557
7377 1989 1 .005867993 27.22292 .001597439
7690 1989 1 .0040071923 -47.93108 -.0019206906
5412 1989 1 .005992175 -26.92 -.0016130934
5709 1989 1 .010813965 38.9516 .0042122123
7643 1989 1 .09222806 5.089854 .004694274
3671 1989 1 .012384918 -22.748404 -.002817371
7206 1989 1 .031833626 28.65433 .009121712
7390 1989 2 .005273794 -12.060838 -.0006360637
7264 1989 2 .005712387 -25.24261 -.0014419556
6205 1989 2 .0021875934 -12.025 -.00026305811
5855 1989 2 .02129368 -33.009 -.00702883
\
end

[Jesse Tielens]

Just to clarify some information, after reading this post again I think I did not make myself very clear (but I can't edit it again).

I currently calculate std. dev and T-score of a yearly portfolio constructed on CI_quintile by doing this:

//calculate std. dev and T-score for zero-mean weighted abnormal return portfolios
gen wx = (weight * ((AR - 0)^2))
bysort year CI_quintile: egen wx_sum = sum(wx)
gen std_dev = sqrt((wx_sum/total_weight))
bysort year CI_quintile: egen portfolio_weighted_mean = mean(weight * AR)
bysort year CI_quintile: gen t_score = (0 - portfolio_weighted_mean)/(sqrt(std_dev / _N))

Which works, but is time-consuming and it seems to me like there should be a shortcut to doing this. If I can find a way to have Stata calculate the T-scores and standard deviations, this also allows me to use clustered errors for example. I just can't figure that shortcut out. For example, if I have Stata do:

gen wAR = weight * ARIf
bysort CI_quintile year: ttest wAR == 0

This produces a very different result than what I calculated myself.

[Carole J. Wilson]

Unless I am missing something, there are a few mistakes in your implementation of the weighted standard deviation formula

1) you assume zero mean within year & quintile, but at least in this sample, that's not true
2) you do not use (m-1)/m in the denominator. In your larger dataset, this might be negligible, but here it isn't

The Stata's results using -aweights- produces standard deviations very close to the corrected code (I get identical weighted means by hand). Perhaps you can use -aweights- to get your t-scores, etc

bysort year CI_quintile: egen double mean=mean(AR)
gen double wx  = (weight * ((AR - mean)^2))
bysort year CI_quintile: egen double wx_sum = sum(wx)
bysort year CI_quintile: egen double total_weight=total(weight)
bysort year CI_quintile: gen groupn=_N
gen groupn_1=groupn-1
gen double n1divn=groupn_1/groupn
gen double std_dev = sqrt((wx_sum/(total_weight*n1divn)))
egen tag=tag(year CI_quintile)
list year CI_quintile std_dev if tag==1

bysort year : tab CI_quintile [aw=weight] ,sum( AR)

[Jesse Tielens]

I made a mistake and this got posted twice. I apologize for the confusion that this may have caused. I suggest we continue all discussion related to my question in this thread. Clyde Shechter provided some good input in the other thread, which we can just leave to sink to the inactive part of the forum I'd say.

Clyde Schechter said: Without having delved too deeply into it, I think the main problem is that you have used the -semean- statistic in your -tabstat- command. -semean- calculates the standard error of the mean, which is not the standard deviation. To get the standard deviation, use -sd- in your -statistics()- option, not -semean-.

Also, it may not be necessary to dance around the weighting by explicitly calculating wAR and descriptive statistics from there. The -tabstat- command accepts -aweights-, which may give you what you are looking for. If not, look, instead at using -summarize-, which accepts -iweights-. One of these is likely to be what you want, with less work on your part to get it.

Thank you both for your help in this matter.

1) you assume zero mean within year & quintile, but at least in this sample, that's not true
2) you do not use (m-1)/m in the denominator. In your larger dataset, this might be negligible, but here it isn't

I think both of these problems are solved in the entire sample. I just checked and mean(AR) for the entire population is 0.0003.

Both of you suggest using aweight, which I think might very well be the solution. The edited code provided by Carole J. Wilson does indeed produce the same standard deviations as tabstat with aweights. I'm having trouble understanding how aweight works though, for example I can do:

tabstat AR [aw=weight], by(CI_quintile) statistics(sum count mean sd median max min) columns(statistics)

Producing the following result:

Does this imply that each of the stats is multiplied by the observations' weight before calculating it? So, I that'd mean for example that the column 'sum' is the sum of (weight * AR)?

Or should i do 'tabstat wAR', so as to tell Stata that each of the observations was weighted by 'weight'?

[Steve Samuels]

I want to point out a discrepancy between your sample data and your weight formula. In the first observation: weight, AR, and wAR are, respectively

.003308191 .04755842 1.5733234e-06

The product of weight & AR is

003308191*.04755842 = 0.00015733234 = 1.5733234 e-04

Something is off by a factor of 100.

[Jesse Tielens]

You are right indeed Steve, very attentive of you!

On a later line in the DO-file I had multiplied AR by 100, so it'd display nicely as a percentage in the outputted table. But this does cause problems with the weights, etc. later on.

wAR = 1.577 e-06

Was correct, but AR was originally 0.004756.

For those interested, I have included a correct dataex here.

[Steve Samuels]

I know nothing about financial statistics, but I question this whole project as it is now described.

1. With an average of N firms per year in each quintile, you will be doing 5N tests per year. With N = 40 (say) and 2 years, that would be 400 tests. The likelihood of falsely declaring a return as "significant" will be very high. You need a multiple comparisons adjustment to those tests

2. The t-test was designed to test hypotheses about means. You are testing hypotheses about single observations of unknown distributions. I strongly doubt that a Gaussian or t-distribution is the proper reference distribution. Therefore, even the single-test Type I error rates are likely to be far from a nominal value like 0.05.

3. It's true that t-tests were once used to detect outliers (assuming a normal distribution). However extreme values inflate both means and standard deviations. The t-statistics are shrunk, reducing the likelihood that those extreme values will be flagged. Such t- or Guassian-based methods have been discredited and robust methods are now preferred

5. If there will be costs to making bad predictions and benefits to making good ones, you are probably better off with a (robust) decision procedure that takes those costs and benefits into account.

[Jesse Tielens]

I can understand your concerns, as a post on Statalist is in no way detailed enough to explain my entire methodology. But I assure you there's an extensive methodology behind all this, which is tried and tested by well-known scholars. I haven't thought of this methodology myself, as I'm still young. It's hard to explain what you're doing over the internet though, you also want to keep it short in order to get some help.

I'm just reading books and papers and using their methodology which is vetted by other researchers. For example the methodology I'm following is specified in 'Capital investment and stock returns' by S. Titman and cited well over 3,000 times. But, they just report T-scores and assume every knows how they've obtained those, I'm unsure how to get those for this specific type of dataset.

1. With an average of N firms per year in each quintile, you will be doing 5N tests per year. With N = 40 (say) and 2 years, that would be 400 tests. The likelihood of falsely declaring a return as "significant" will be very high. You need a multiple comparisons adjustment to those tests

On the contrary, each year has 5 portfolio's compromised of thousands of firms. What I'm interested in is whether the value-weighted mean return of each portfolio is different. I have data on 29 years, with 5 portfolios. That makes an N of only 125 in the end. The whole thing I don't get about this, is that if you're going to do a T-test on observations that are means in itself. Shouldn't the variance of those observations you used to calculate the mean with be preserved.

E.g. If I suspect that the firms in CI_quintile 1 have higher AR, I look at the data and see that the mean of CI_quintile 1 in 1999 is high, but has a very high standard deviation as well. I don't think it's fair to use that group's mean as an observation in the T-test of means. Because that group's mean isn't a fair reflection of the underlying group's data. That is why I want to investigate the std. dev and T-scores of the yearly means.

The t-test was designed to test hypotheses about means. You are testing hypotheses about single observations of unknown distributions. I strongly doubt that a Gaussian or t-distribution is the proper reference distribution. Therefore, even the single-test Type I error rates are likely to be far from a nominal value like 0.05.

AR as a whole has a mean of 0.003, so very close to zero and its histogram looks fairly normally distributed actually in my opinion. With the exception that there are a few extreme outliers at the right. But, most methodology papers I've read say that those outliers aren't a problem for sample sizes this large.


Known problems are that the data is generally a bit heteroskedastic (as is often the case with finance) and skewed to the right because of some extreme outliers at the end (also often the case in finance). Now, I can do some RObustness testing for that. Most often people do a Kruskal - Wallis test between groups as Robustness, because that can be done on any distribution. Essentially testing to see if their conclusion still holds up, even if the distribution isn't assumed to be normal.

It's true that t-tests were once used to detect outliers (assuming a normal distribution). However extreme values inflate both means and standard deviations. The t-statistics are shrunk, reducing the likelihood that those extreme values will be flagged. Such t- or Guassian-based methods have been discredited and robust methods are now preferred

Wouldn't the large sample size overcome these problems?

[Steve Samuels]

I'm unclear about the ultimate analysis. Suppose that you reject none of the means. What's next?

[Steve Samuels]

I understand now that the interest is in quintile means. I'm obviously approaching this problem in ignorance of standard methods in the area.

1. Is there reason to think that quintile meansm and standard deviations should be similar? If not, I don't see the logic in omitting a quintile because it looks different.

2. Comparisons of one quintile against the whole are actually comparisons against the other quintiles. See https://stats.stackexchange.com/questions/52614/are-there-any-statistically-significant-differences-between-my-study-population/ The sizes of t-scores are reduced because each quintile contributes to the whole. As a technical matter, the t-scores do not reflect variation in the second group, thus are too big. It appears that a two-sample t-score is called for. Note that the (two-sample?) KW test , while insensitive to outliers, 1) does not compare means; and 2) is not robust to differences in standard deviations.

3. If interest is in means, I suggest a trimmed-mean as a fair reflection of the vast majority of the underlying group's data.

4. To answer the last question: if the fraction of extreme points is pretty constant, than large sample size doesn't solve the problem.

[Steve Samuels]

I should have added that a two-sample t-test or KW test comparing one quintile to the whole would be invalid. The assumption of both tests is that the two groups are independent, not true of the part-whole groups.