Last data update: 2014.03.03

R: Calculate a credible interval from a numerically specified...
credIntNumR Documentation

Calculate a credible interval from a numerically specified posterior CDF

Description

Calculates a lower, upper, or two-sided credible interval from the numerical posterior CDF.

Usage

credIntNum(theta, cdf, conf = 0.95, type="twosided")

Arguments

theta

the values over which the the posterior CDF is specified

cdf

the values of the CDF, F(θ) = int_{-∞}^{θ}f(t).df where f(t) is the PDF.

conf

the desired 'confidence' level

type

the type of interval to return, 'lower' = one sided lower bound, 'two-sided' = two - sided, or 'upper' = one sided upper bound. It is sufficient to use 'l','t' or 'u'

Details

This function uses linear interpolation to calculate bounds for points that may not be specified by CDF

Value

a list containing the elements lower.bound, uppper.bound or both depending on type

Examples

## commands for calculating a numerical posterior CDF.
## In this example, the likelihood is proportional to
## eqn{\theta^{3/2}\times exp(-\theta/4)} and a N(6, 9) prior is used.
theta <- seq(from = 0.001, to = 40, by = 0.001)
prior <- dnorm(theta,6,3)
ppnLike <- theta^1.5*exp(-theta/4)
ppnPost <- prior*ppnLike
scaleFactor <- sintegral(theta, ppnPost)$int
posterior <- ppnPost/scaleFactor
cdf <- sintegral(theta, posterior)$y
ci<-credIntNum(theta, cdf)
par(mfrow=c(2,2))
plot(prior ~ theta, type = 'l',  main = "Prior N(6, 9)")
plot(ppnLike ~ theta, type = 'l', main = "Proportional likelihood")
plot(posterior ~ theta, type = 'l', main = "Posterior")
abline(v=c(unlist(ci)))

Results


R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(Bolstad2)
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/Bolstad2/credIntNum.Rd_%03d_medium.png", width=480, height=480)
> ### Name: credIntNum
> ### Title: Calculate a credible interval from a numerically specified
> ###   posterior CDF
> ### Aliases: credIntNum
> 
> ### ** Examples
> 
> ## commands for calculating a numerical posterior CDF.
> ## In this example, the likelihood is proportional to
> ## eqn{\theta^{3/2}\times exp(-\theta/4)} and a N(6, 9) prior is used.
> theta <- seq(from = 0.001, to = 40, by = 0.001)
> prior <- dnorm(theta,6,3)
> ppnLike <- theta^1.5*exp(-theta/4)
> ppnPost <- prior*ppnLike
> scaleFactor <- sintegral(theta, ppnPost)$int
> posterior <- ppnPost/scaleFactor
> cdf <- sintegral(theta, posterior)$y
> ci<-credIntNum(theta, cdf)
Credible interval is : (2.02114790458962,11.4246333371552)
> par(mfrow=c(2,2))
> plot(prior ~ theta, type = 'l',  main = "Prior N(6, 9)")
> plot(ppnLike ~ theta, type = 'l', main = "Proportional likelihood")
> plot(posterior ~ theta, type = 'l', main = "Posterior")
> abline(v=c(unlist(ci)))
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>