0
$\begingroup$

The documentation for PDF has an example for plotting a confidence interval of a standard normal distribution:

pdf = PDF[NormalDistribution[], x]
Show[Plot[pdf, {x, -1.7, -1}, Filling -> Axis, FillingStyle -> Green],
  Plot[pdf, {x, -1, 1}] , 
 Plot[pdf, {x, 1, 1.7}, Filling -> Axis, FillingStyle -> Green], 
 PlotRange -> {0, 0.4}, AxesOrigin -> {0, 0.1}, 
 Ticks -> {Automatic, None}]

The plot range is set such that this code would not work in a generalized case for NormalDistribution[u,s].

What is the appropriate format for {x, xmin, xmax} in the Plot function above for the general case of a normal distribution with mean = u and standard deviation = s?

$\endgroup$
2
  • $\begingroup$ @b.gatessucks That's what I'm doing, my question is how to programmatically determine the correct xmin/xmax for any given u and s. It's more of a stats question than a Mathematica question and may not be appropriate here. $\endgroup$ Nov 15, 2013 at 19:16
  • $\begingroup$ Is it bad etiquette to write notes to myself in the comment section? I know when I'm looking for this Q&A I'll undoubtedly start by looking at this question. $\endgroup$ Nov 15, 2013 at 20:26

3 Answers 3

2
$\begingroup$

This code does it

pdf = PDF[NormalDistribution[2, 1], x]
With[{mu = 2, s = 1},
 Show[
  Plot[pdf, {x, mu - 2 s, mu - s}, Filling -> Axis, 
   FillingStyle -> Green],
  Plot[pdf, {x, mu - s, mu + s}],
  Plot[pdf, {x, mu + s, mu + 2 s}, Filling -> Axis, 
   FillingStyle -> Green], PlotRange -> {0, 0.4}, 
  AxesOrigin -> {0, 0.1}, Ticks -> {Automatic, None}
  ]
 ]

Essentially, you are plotting 3 intervals $x \in [\mu - \sigma, \mu + \sigma]$, $x \in [\mu - 2 \sigma, \mu - \sigma]$, and $x \in [\mu + \sigma, \mu + 2 \sigma]$ where the last two intervals are also shaded green.

$\endgroup$
1
  • $\begingroup$ Thanks for the answer and the stats lesson. $\endgroup$ Nov 15, 2013 at 19:50
1
$\begingroup$

The limits of a 100c% confidence interval for a normal variable with mean u and standard deviation s are u ± s*InverseErf[c]*Sqrt[2]. The example doesn't say how the 1.7 cutoff was chosen, but Erf[1.7/Sqrt[2]] = .910869, so the entire plot is only a 91% confidence interval for a standard normal, not the almost-100% that is usual. To give a correct visual impression, the green should have been extended to ±3 or so with the x-axis at y=0. The white area, whose limits area ±1, is a 68% interval for a standard normal. This is how they should have done it to show a 68% interval:

Show[
  Plot[pdf, {x, -3.3, -1}, Filling -> Axis, FillingStyle -> Green],
  Plot[pdf, {x, -1, 1}, Filling -> None],
  Plot[pdf, {x, 1, 3.3}, Filling -> Axis, FillingStyle -> Green],
  PlotRange -> {0, .4}, Axes -> {True, False}, Ticks -> {Automatic, None}] 

enter image description here

$\endgroup$
1
  • $\begingroup$ $0.7$ is approx. $0.68$, so they're using one standard deviation. $\endgroup$
    – rcollyer
    Nov 15, 2013 at 21:06
0
$\begingroup$

I am not sure whether this will be useful. If the intention is to display the standard symmetric confidence interval of an arbitrary normal distribution, then the following will do this. Noting:

  1. the confidence interval is 'highlighted' by the complement of the filled regions. If the intention is to highlight the interval in a special colour just reverse the RegionFunction inequality.
  2. I have used Filling to 0. as Filling to Axis has some issues which I do not fully understand.

ciplot[m_, s_, ci_, opts : OptionsPattern[]] :=
  Module[{pdf, q, ns, p1, p2},
   pdf = NormalDistribution[m, s];
   q = Quantile[pdf, (1 + ci)/2];
   ns = (q - m)/s;
   p1 = Plot[PDF[pdf, x], {x, m - 1.5 ns s, m + 1.5 ns s}]; 
   p2 = Plot[PDF[pdf, x], {x, m - 1.5 ns s, m + 1.5 ns s}, 
     Filling -> 0., Evaluate@FilterRules[{opts}, Options[Plot]], 
     RegionFunction -> Function[x, Abs[x - m] > ns s]];
   Show[p2, p1, AxesOrigin -> {0, 0}, PlotRange -> All]];

I generate the animated gif for illustration purposes and manipulate could be used varying mean and standard deviation.

enter image description here

$\endgroup$

Your Answer

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

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