Fitting standard deviation of a normal distribution from one or more quantiles and known mean

Description

get.norm.sd returns the standard deviation of a normal distribution where the pth percentiles match with the quantiles q.

Usage

1
2
get.norm.sd(p = c(0.025, 0.5, 0.975), q, show.output = TRUE, plot = TRUE,
                   fit.weights = rep(1, length(p)), scaleX = c(0.1, 0.9), ...)

Arguments

p

numeric, vector of probabilities.

q

numeric, vector of quantiles corresponding to p.

show.output

logical, if TRUE the optim result will be printed (default vaule is TRUE).

plot

logical, if TRUE the graphical diagnostics will be plotted (default value is TRUE).

fit.weights

numerical vector of the same length as a probabilities vector p containing positive values for weighting quantiles. By default all quantiles will be weighted by 1.

scaleX

numerical vector of the length 2 containing values (from the open interval (0, 1)) for scaling quantile-axis (relevant only if plot = TRUE). The smaller the left value, the further the graph is extrapolated within the lower percentile, the greater the right value, the further it goes within the upper percentile.

...

further arguments passed to the functions plot and points (relevant only if plot = TRUE).

Details

The number of probabilities and the number of quantiles must be identical and should be at least two. get.norm.sd uses the central limit theorem and the linear regression.

If show.output = TRUE the output of the function lm will be shown.

The items of the probability vector p should lie between 0 and 1.

The items of the weighting vector fit.weights should be positive values.

The function will be meaningful only if the quantile comes from a normal distribution.

Value

Returns an estimated standard deviation or missing value

Note

It should be noted that the data must be normally distributed, or the central limt theorem must hold for large (enough) samples sizes.

Author(s)

Matthias Greiner matthias.greiner@bfr.bund.de (BfR),
Katharina Schueller schueller@stat-up.de (STAT-UP Statistical Consulting),
Natalia Belgorodski belgorodski@stat-up.de (STAT-UP Statistical Consulting)

See Also

See pnorm for distribution implementation details.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
q <- stats::qnorm(p = c(0.025, 0.5, 0.975), mean = 0, sd = 2)
old.par <- graphics::par(mfrow = c(2, 3))
get.norm.sd(q = q)
get.norm.sd(q = q, scaleX = c(0.0001, 0.9999))
get.norm.sd(q = q, fit.weights = c(10, 1, 10))
get.norm.sd(q = q, fit.weights = c(1, 10, 1))
get.norm.sd(q = q, fit.weights = c(100, 1, 100))
get.norm.sd(q = q, fit.weights = c(1, 100, 1))
graphics::par(old.par)

q <- stats::qnorm(p = c(0.025, 0.5, 0.975), mean = 176, sd = 15)
old.par <- graphics::par(mfrow = c(2, 3))
get.norm.sd(q = q)
get.norm.sd(q = q, fit.weights = c(10, 1, 10))
get.norm.sd(q = q, fit.weights = c(1, 10, 1))
get.norm.sd(q = q, fit.weights = c(100, 1, 100))
get.norm.sd(q = q, fit.weights = c(1, 100, 1))
graphics::par(old.par)

## The estimation model is not suitable for the following quantiles.
## Because the quantile is unsymmetrical, which could not be from a normally distributed data.
q <- c(-2, 30, 31)
old.par <- graphics::par(mfrow = c(2, 3))
get.norm.sd(q = q)
get.norm.sd(q = q, fit.weights = c(10, 1, 10))
get.norm.sd(q = q, fit.weights = c(1, 10, 1), scaleX = c(0.0001, 0.9999))
get.norm.sd(q = q, fit.weights = c(100, 1, 100))
get.norm.sd(q = q, fit.weights = c(1, 100, 1), scaleX = c(0.0001, 0.9999))
graphics::par(old.par)

## Estimating from actually exponentially distributed data
x.exp <- rexp(n = 10, rate = 5)
mean(x.exp)
stats::sd(x.exp)
q <- quantile(x.exp, c(0.025, 0.5, 0.975))
old.par <- graphics::par(mfrow = c(2, 3))
get.norm.sd(q = q)
get.norm.sd(q = q, fit.weights = c(1, 10, 1))
get.norm.sd(q = q, fit.weights = c(10, 1, 10))
get.norm.sd(q = q, fit.weights = c(1, 100, 1))
get.norm.sd(q = q, fit.weights = c(100, 1, 100))
graphics::par(old.par)

## other examples
q <- stats::qnorm(p = c(0.025, 0.5, 0.975), mean = 1, sd = 1)
get.norm.sd(q = q)

q <- stats::qnorm(p = c(0.025, 0.5, 0.975), mean = 1, sd = 0.5)
get.norm.sd(q = q)

q <- stats::qnorm(p = c(0.025, 0.5, 0.975), mean = 0.01, sd = 0.1)
get.norm.sd(q = q)

Want to suggest features or report bugs for rdrr.io? Use the GitHub issue tracker.