# NOT RUN {
## Example 1 - product of 2 probability parameters for which only the 95% CIs are reported
dist1<-getBetaFromCI(qLow=0.4,qUpp=0.6,alpha=0.05)
dist2<-getBetaFromCI(qLow=0.7,qUpp=0.9,alpha=0.05)
distListEx<-list(dist1$r,dist2$r)
combFunEx<-function(pars){pars[[1]]*pars[[2]]}
bootComb(distList=distListEx,
combFun=combFunEx,
doPlot=TRUE,
method="hdi",
N=1e5, # reduced from N=1e6 so that it runs quicker; larger values => more accurate
seed=352)
# Alternatively, the same example can be run in just 2 lines of code:
combFunEx<-function(pars){pars[[1]]*pars[[2]]}
bootComb(distributions=c("beta","beta"),
qLowVect=c(0.4,0.7),
qUppVect=c(0.6,0.9),
combFun=combFunEx,
doPlot=TRUE,
method="hdi",
N=1e5, # reduced from N=1e6 so that it runs quicker; larger values => more accurate
seed=352)
## Example 2 - sum of 3 Gaussian distributions
dist1<-function(n){rnorm(n,mean=5,sd=3)}
dist2<-function(n){rnorm(n,mean=2,sd=2)}
dist3<-function(n){rnorm(n,mean=1,sd=0.5)}
distListEx<-list(dist1,dist2,dist3)
combFunEx<-function(pars){pars[[1]]+pars[[2]]+pars[[3]]}
bootComb(distList=distListEx,combFun=combFunEx,doPlot=TRUE,method="quantile")
# Compare with theoretical result:
exactCI<-qnorm(c(0.025,0.975),mean=5+2+1,sd=sqrt(3^2+2^2+0.5^2))
print(exactCI)
x<-seq(-10,30,length=1e3)
y<-dnorm(x,mean=5+2+1,sd=sqrt(3^2+2^2+0.5^2))
lines(x,y,col="red")
abline(v=exactCI[1],col="red",lty=3)
abline(v=exactCI[2],col="red",lty=3)
## Example 3 - same as Example 1 but assuming the 2 parameters to be dependent / correlated
combFunEx<-function(pars){pars[[1]]*pars[[2]]}
bootComb(distributions=c("beta","beta"),
qLowVect=c(0.4,0.7),
qUppVect=c(0.6,0.9),
Sigma=matrix(byrow=TRUE,ncol=2,c(1,0.5,0.5,1)),
combFun=combFunEx,
doPlot=TRUE,
method="hdi",
N=1e5, # reduced from N=1e6 so that it runs quicker; larger values => more accurate
seed=352)
# }
Run the code above in your browser using DataLab