Nothing
## MRRln2p.r file
##
## Author: Jacob T. Ormerod
## (c)2014 OpenReliability.org
##
MRRln2p<-function(x, s=NULL, bounds=FALSE, CI=0.90, show=FALSE) {
xName<-paste(deparse(substitute(x),500),collapse = "\n")
thisTitle=paste(xName,"lognormal fit")
gotppp<-getPPP(x, s)
probability<-as.vector(gotppp[,2])
x<-as.vector(gotppp[,1])
fit<-lslr(gotppp, dist="lnorm")
fit<-c(fit,LL=LLln(x,s,fit[1],fit[2]))
if(bounds==TRUE) {
## get the parameters to initiate final pivotal analysis
## using an adaptation of Jurgen Symynck's method for weibull applied to lognormal
completePPP<-getPPP(x)
comp_prob<- as.vector(completePPP[,2])
P1=mean(qnorm(probability,0,1))
P2=sd(qnorm(probability,0,1))/sd(qnorm(comp_prob,0,1))
## division by sd(qnorm(comp_prob,0,1))
## corrects for resolution of sd for complete failures to 1.0
## critically needed for lognormal pivotal parameter determination
## descriptive quantiles for comparison with SuperSMITH (limit of 15 values)
dq<-c(.01, .02, .05, .10, .15, .20, .30, .40, .50, .60, .70, .80, .90, .95, .99)
pivotals<-pivotalMC(gotppp,dist="lnorm",R2=0, CI=CI,unrel=dq,P1=P1,P2=P2)
## check the slope of the median pivotals to get correction to 1.0
median_slope<-(qnorm(dq[15], 0, 1) - qnorm(dq[1], 0, 1))/(pivotals[15,2]-pivotals[1,2])
median_intercept<-pivotals[10,2]-qnorm(dq[10], 0, 1)/median_slope
adj_piv<-(pivotals-median_intercept)*median_slope
## intermediate readings of the slope were generated during development, now not used
## median_slope2<-(qnorm(dq[15], 0, 1)-qnorm(dq[9], 0, 1))/(adj_piv[15,2]-adj_piv[9,2])
## median_intercept2<-adj_piv[10,2]-qnorm(dq[10], 0, 1)/median_slope2
## interpret the pivotals for the log plot
plot_piv<-(adj_piv)*fit[2]+fit[1]
## again intermediate readings of the slope and intercept were generated only for development
## median_slope3<-(qnorm(dq[15], 0, 1)-qnorm(dq[9], 0, 1))/(plot_piv[15,2]-plot_piv[9,2])
## median_intercept3<-plot_piv[10,2]-log(log(1/(1-dq[10])))/median_slope3
## prepare the pivotals for plotted output
LB<-as.vector(plot_piv[,1])
DATUM<-as.vector(plot_piv[,2])
HB<-as.vector(plot_piv[,3])
## prepare the pivotals for print output
print_piv<-exp(plot_piv)
outDF<-data.frame(unrel=dq,
LB=as.vector(print_piv[,1]),
DATUM=as.vector(print_piv[,2]),
HB=as.vector(print_piv[,3]))
}
if(show==TRUE) {
plot(log(x),qnorm(probability),pch=19,col="red",
main=thisTitle)
Xintercept<--fit[1]/fit[2]
abline(Xintercept,1/fit[2],col="blue")
if(bounds==TRUE) {
lines(LB,qnorm(dq),col="blue")
lines(HB,qnorm(dq),col="blue")
lines(DATUM, qnorm(dq), col="magenta")
}
}
if(bounds==TRUE) {
return(list(fit,outDF))
}else{
return(fit)
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.