9
$\begingroup$

I have a dataset for which Mathematica easily creates a Histogram. However I also need Mathematica to show error bars corresponding to 3 standard deviations for each bar. This is something similar to the ErrorBar function for BarChart. I cannot find any way to do this.

The difference between my case and the ErrorBar package avaliable (http://reference.wolfram.com/language/ErrorBarPlots/ref/ErrorBar.html) is that it requires data for the magnitude of error for each bin/bar. In my case, I am giving Mathematica a raw data file to plot into a Histogram. It should be able to calculate the error by itself. There is no way in which I can give Mathematica the bar heights, bin boundaries and corresponding errors from which it will give me a histogram. If someone can point that out, even that would be a potential solution to my problem. The ErrorBar method works for BarCharts, not Histograms as I understand. Could anyone help?

Update:

Using the code @kglr helped me with I could implement the error bars in my Histogram. It is a Histogram that takes its y-values as {"Log","Count"}. So some modification was required.

 ceF[d_: .2, nsd_: 3, color_: Automatic][cedf_: "Rectangle"] := 
     Module[{e = 
         nsd Sqrt[Exp[#[[2, 2]]]]}, {ChartElementData[cedf][##], Thick, 
        color /. 
         Automatic -> Darker[Charting`ChartStyleInformation["Color"]],
        If[Exp[#[[2, 2]]] - e != 0,
         {Line[{{Mean@#[[1]], Log[Exp[#[[2, 2]]] - e]}, {Mean@#[[1]], 
             Log[Exp[#[[2, 2]]] + e]}}], 
          Line[{{#[[1, 1]] + d/2, 
             Log[Exp[#[[2, 2]]] - e]}, {#[[1, 2]] - d/2, 
             Log[Exp[#[[2, 2]]] - e]}}], 
          Line[{{#[[1, 1]] + d/2, 
             Log[Exp[#[[2, 2]]] + e]}, {#[[1, 2]] - d/2, 
             Log[Exp[#[[2, 2]]] + e]}}]},
         {Line[{{Mean@#[[1]], 0}, {Mean@#[[1]], 
             Log[Exp[#[[2, 2]]] + e]}}], 
          Line[{{#[[1, 1]] + d/2, 0}, {#[[1, 2]] - d/2, 0}}], 
          Line[{{#[[1, 1]] + d/2, 
             Log[Exp[#[[2, 2]]] + e]}, {#[[1, 2]] - d/2, 
             Log[Exp[#[[2, 2]]] + e]}}]}
         ]
        }] &

Incidentally the highest bin of my histogram faces this issue where Exp[#[[2, 2]]] - e is exactly zero. This modification should have been able to solve the issue but what I get is this, enter image description here

and I don't know why the Bars dont start from the zero of the axis. That is to say, the bars should be starting from the y-tic marked "1", but they are starting from somewhere above "3." Why is that so? What do I do to make it start from "1"?

Incidentally, for the people who think this is a duplicate of the linked Q/A's, I'm sorry I could not make the connection since I am not an expert at Mathematica. However, the last time I looked it was not called Mathematica StackExchange for Experts Only. Please have this much respect for another person to realise that he would not have posted the question explaining why it is not a duplicate if he could solve the problem from there. I am a professional and my time and energy has some value. I would not be here if I could solve the problem by myself.

$\endgroup$
7
  • 1
    $\begingroup$ When bar height is "Count" (or "PDF") , how would you calculate "3 standard deviations for each bar"? $\endgroup$
    – kglr
    Jul 15, 2018 at 13:41
  • $\begingroup$ When the count is N the standard deviation is Sqrt[N]. You can go through this link for an explanation. suchideas.com/articles/maths/applied/histogram-errors $\endgroup$ Jul 16, 2018 at 6:38
  • 1
    $\begingroup$ rhermans, Jens, Henrik, MarcoB, @halirutan, imho, this question is not a duplicate of the linked Q/A. The OP is quite clear that he wants a Histogram, and none of the answers in the linked q/a works without additional processing. $\endgroup$
    – kglr
    Jul 16, 2018 at 23:37
  • 1
    $\begingroup$ @kglr I looked over it and decided that the linked QA is sufficient to come up with a solution and since you already gave an answer, that should be enough help for the user. I'm sure I still judge questions differently in this regard though and my decision was definitely supported by the existing 4 close votes. Nevertheless, if others feel as you do and there is another open vote, I'm all in for reopening this. $\endgroup$
    – halirutan
    Jul 16, 2018 at 23:45
  • 1
    $\begingroup$ @HeerakBanerjee "why the Bars don't start from the zero of the axis" You may need a slightly bigger monitor given how far down is zero in Log scale.... :) $\endgroup$
    – rhermans
    Jul 17, 2018 at 8:58

2 Answers 2

16
$\begingroup$

For the case where the height function is "Count", we can use the formula from the linked page in a custom ChartElementFunction with the sample size (Length[data]) passed as metadata:

ceF[d_: .2, nsd_: 3, color_: Automatic][cedf_: "Rectangle"] := 
  Module[{e = nsd /2 Sqrt[#[[2, 2]] (1 - #[[2, 2]]/ #3[[1]])]}, 
   {ChartElementData[cedf][##], Thick, 
     color /. Automatic -> Darker[Charting`ChartStyleInformation["Color"]], 
     Line[{{Mean@#[[1]], #[[2, 2]] - e}, {Mean@#[[1]], #[[2, 2]] + e}}], 
     Line[{{#[[1, 1]] + d/2, #[[2, 2]] - e}, {#[[1, 2]] - d/2, #[[2, 2]] - e}}], 
     Line[{{#[[1, 1]] + d/2, #[[2, 2]] + e}, {#[[1, 2]] - d/2, #[[2, 2]] + e}}]}]&

Examples:

SeedRandom[1]
data = RandomVariate[NormalDistribution[0, 1], 200];
Histogram[data -> Length[data], ChartStyle -> 43, ChartElementFunction -> ceF[][]]

enter image description here

Histogram[data -> Length@data, ChartStyle -> 43, 
 ChartElementFunction -> ceF[.2, 3, Black]["GlassRectangle"]]

enter image description here

Update: to make it work with non-default BarOrigin settings:

ClearAll[ceF]
ceF[d_: .2, nsd_: 3, color_: Automatic][cedf_: "Rectangle"] := 
 Module[{bo = Charting`ChartStyleInformation["BarOrigin"], 
     col = Darker[Charting`ChartStyleInformation["Color"]], box = #, tf, e},
    tf = Switch[bo, Left | Right, Reverse, _, Identity];
    box = Switch[bo, Bottom, box, Top, {box[[1]], Reverse[box[[2]]]}, Left,
       Reverse@box, Right, {box[[2]], Reverse@box[[1]]}];
    e = nsd /2 Sqrt[Abs@box[[2, 2]] (1 - Abs@box[[2, 2]]/#3[[1]])]; 
  {ChartElementData[cedf][##], Thick, color /. Automatic -> col, 
    Line[tf /@ {{Mean@box[[1]], box[[2, 2]] - e},
       {Mean@box[[1]], box[[2, 2]] + e}}], 
    Line[tf /@ {{box[[1, 1]] + d/2, box[[2, 2]] - e},
       {box[[1, 2]] - d/2, box[[2, 2]] - e}}], 
    Line[tf /@ {{box[[1, 1]] + d/2, box[[2, 2]] + e}, 
       {box[[1, 2]] - d/2, box[[2, 2]] + e}}]}] &

Example:

Grid[Partition[Histogram[data -> Length@data, ChartStyle -> 43, 
     ChartElementFunction -> ceF[][], ImageSize -> 300, 
     BarOrigin -> #] & /@ {Bottom, Top, Left, Right}, 2], 
 Dividers -> All]

enter image description here

$\endgroup$
2
  • $\begingroup$ Maybe my memory is failing me again but I don't recall seeing the syntax data -> Length[data] before; is this documented? $\endgroup$
    – Mr.Wizard
    Jul 17, 2018 at 6:42
  • 1
    $\begingroup$ @Mr.Wizard, $data -> metadata$ is the standard way to attach metadata for *Charts. It is not documented in Histogram pages but it works. Also not documented but wrappers (Button, Hyperlink, ... PopupWindow etc) can be used with Histogram input data. $\endgroup$
    – kglr
    Jul 17, 2018 at 6:55
4
$\begingroup$

This answer below is not directly what you asked but rather about what you should consider doing. With more than 50 or so data points you should consider avoiding histograms completely. More often than not you probably envision some smooth density function that you're trying to estimate. Further, adding in error bars makes for a very messy and maybe even difficult to interpret figure.

And finally taking the log or square root for displaying the count or density makes no sense. That destroys the feature of the area under the histogram (or density curve) "summing" to 1 or summing to the sample size and makes comparisons among different datasets risky at best.

Taking the square root of the count can get you an estimate of the standard error associated with the count - formally called the rootogram. My complaint is about then being losing the ability to make sense of comparing different datasets.

Transforming the data (i.e., the raw data and not the counts is a different and perfectly fine approach.

Using a nonparametric density estimate with bootstrap-created confidence bands might show the features of your data much better. Consider the following:

(* Generate some data *)
SeedRandom[1]
n = 200;
data = RandomVariate[NormalDistribution[0, 1], n];

(* Estimate density function *)
skd = SmoothKernelDistribution[data];

(* Determine some bounds to evaluate the density function *)
sd = StandardDeviation[data];
xmin = Min[data] - sd;
xmax = Max[data] + sd;

(* Generate bootstrap samples and determine density values along a 
   grid between xmin and xmax *)
nboot = 1000;
ngrid = 100;
densityValues = ConstantArray[0, {nboot, ngrid + 1}];
Do[bootData = RandomChoice[data, n];
 skdboot = SmoothKernelDistribution[bootData];
 densityValues[[iboot, All]] = 
  Table[PDF[skdboot, xmin + (xmax - xmin) i/ngrid], {i, 0, ngrid}],
 {iboot, nboot}]

(* Choose some level for the confidence bands and calculate percentiles *)
confLevel = 0.95;
xvalues = Table[xmin + (xmax - xmin) i/ngrid, {i, 0, ngrid}];
lower = Transpose[{xvalues, Quantile[densityValues[[All, #]],(1 - confLevel)/2] & /@ Range[ngrid + 1]}];
upper = Transpose[{xvalues, Quantile[densityValues[[All, #]], 1 - (1 - confLevel)/2] & /@ Range[ngrid + 1]}];

(* Plot results *)
Show[ListPlot[{upper, lower}, Joined -> True, 
  PlotStyle -> {{Blue, Dotted}}],
  Plot[PDF[skd, x], {x, xmin, xmax}, PlotStyle -> Blue]]

Estimate of density and confidence bands

I think the above figure is much cleaner and informative than having a lumpy histogram with error bars sticking all over the place. We all have computers now. There's no need to do things (like histograms) that were all one could do when computational power was low.

$\endgroup$
6
  • 1
    $\begingroup$ +1 for non-parametric estimates. Histograms are a relic. $\endgroup$
    – rhermans
    Jul 17, 2018 at 17:20
  • 2
    $\begingroup$ @rhermans My wife says something similar about my wardrobe. $\endgroup$
    – JimB
    Jul 17, 2018 at 17:21
  • $\begingroup$ @JimB Thanks a lot for your input. This is indeed a wonderful way for bypassing a histogram. However my data represents values detected by a detector. Say the number of particles carrying energy E (x-values) . And as any other detector it has a minimum sensitivity, which becomes my binning size. So in my case, fitting the data to a continuous density function over the entire x-range would be physically incorrect. Hence the Histogram. I need the bins. And then I would be able to compare this dataset to another dataset for EACH INDIVIDUAL bin. Because that is what I would get from the detector. $\endgroup$ Jul 18, 2018 at 6:38
  • $\begingroup$ @HeerakBanerjee "fitting the data to a continuous density function over the entire x-range would be phyically incorrect." Yes, but...the continuous normal distribution can (under some circumstances) provide a good approximation to a very discrete and non-negative Poisson distribution and few have heartburn over that. Let go of histograms! Let your nonparametric density estimate "describe" and summarize your data as it (and histograms) are not meant to "explain" the data. $\endgroup$
    – JimB
    Jul 18, 2018 at 14:09
  • $\begingroup$ Well @JimB I would happily agree with you if I were trying to fit data that already exists. I am giving projections to a machine here. And the only way IT can take data is in bins. I am restricted by the physics of the issue. However I totally agree with you about whatever you said. $\endgroup$ Jul 20, 2018 at 6:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.