SFRATReport.Rmd

---
title: 'Software Failure and Reliability Assessment Tool: Report'
author: "Author Name"
date: '`r format(Sys.time(), "%Y-%m-%d_%H:%M")`'
output:
  #html_document: default
  pdf_document: default
---
  
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
#opts_chunk$set(tidy.opts=list(width.cutoff=100),tidy=TRUE)
```

```{r, echo=FALSE}
source('./SFRAT/utility/data/Data_Tools.R')##DATA PREPROCESSING
d <- dataset #Input excel file with a single sheet for now
cnames <- colnames(d) # Read column names in the input excel file

tryCatch({                                  #Data conversion depending on the type of the input data
  if("FN" %in% cnames && "IF" %in% cnames && "FT" %in% cnames) {
    FT <- d$FT
    IF <- d$IF
    FN <- d$FN
  } else if("FN" %in% cnames && "IF" %in% cnames) {
    FT <- IF_to_FT(d$IF)
    IF <- d$IF
    FN <- d$FN
  } else if("FN" %in% cnames && "FT" %in% cnames) {
    IF <- FT_to_IF(d$FT)
    FT <- d$FT
    FN <- d$FN
  } else if("T" %in% cnames && "FC" %in% cnames && "CFC" %in% cnames) {
    FC <- d$FC
    CFC <- d$CFC
    FT <- FC_to_FT(d$T,d$FC)
    IF <- FT_to_IF(FT)
    FN <- 1:length(FT)
  } else if("T" %in% cnames && "FC" %in% cnames) {
    CFC <- FC_to_CFC(d$FC)
    FT <- FC_to_FT(d$T,d$FC)
    IF <- FT_to_IF(FT)
    FN <- 1:length(FT)
  } else if("T" %in% cnames && "CFC" %in% cnames) {
    FC <- CFC_to_FC(d$CFC)
    FT <- FC_to_FT(d$T,d$FC)
    IF <- FT_to_IF(FT)
    FN <- 1:length(FT)
  } else{
    print("Upload your input data/file formatted according to the SFRAT guidelines")
  }
}, error = function(error_condition){
  print("Unable to load in data.")
})
d <- data.frame("FN"=FN,"IF"=IF,"FT"=FT)
```

# Tab 1: Select, Apply, and Analyze Data

## Sample of the updated data ('`r SheetName`') in different formats:
`r if (x==1){paste0("The table below shows the first ten points of the input data ",SheetName,". 'FC', 'CFC', 'FT', 'IF', and 'FN' indicates failure counts, cumulative failure counts, failure times, interfailure times, and number of failures respectively.")}`

```{r, echo=FALSE}
if("T" %in% cnames && "FC" %in% cnames && "CFC" %in% cnames) {
  list(data.frame("FT"=head(FT),"IF"=head(IF),"FN"=head(FN)),data.frame("FC"=head(FC), "CFC"=head(CFC)))
}else{
  kable(d[1:10,], caption="First ten points of the input data")
}
```

\newpage
## Cumulative failures
`r if (x==1){paste0("The following figure shows the ",SheetName," data  as the cumulative number of failures (FN) detected as a function of cumulative test time (FT). An increasing trend indicates periods where more faults were detected. Ideally, the cumulative number of failures should level off to a horizontal line, indicating that no new faults have been detected.")}`

```{r, echo=FALSE}
plot(
  FT, FN, type="s",
  xlab="Cumulative test time", ylab="Cumulative number of failures",
  main =  bquote("Cumuative Failures vs. cumulative test time:" ~.(SheetName))
)
```

\newpage
## Times between failures/Interfailure times

`r if (x==1) {paste0("The following figure shows the ",SheetName," times between failures (IF) as a function of cumulative test time (FT). An increasing trend indicates periods where fewer faults were detected. Ideally, the time between failures should increase, indicating that no new faults have been detected.")}`

```{r, echo=FALSE}
plot(FT, IF, type="s", xlab="Cumulative test time", ylab="Times between failures",main = bquote("Interfailure times vs. cumulative test time:" ~.(SheetName)))
```

\newpage
## Failure intensity
`r if (x==1){paste0("The following figure shows the ",SheetName," data as the number of failures detected per unit time as a function of cumulative test time (FT). A decreasing trend indicates periods where fewer faults were detected. Ideally, the failure intensity should decrease, indicating that no new faults have been detected.\n\n
A decrease in failure intensity indicates increase in reliability of the software subjected to testing. Ideally, the failure intensity should go to zero.")}`

```{r, echo=FALSE}
plot(FT, 1/IF, type="s",xlab="Cumulative test time",ylab="Number of failures per unit time",main = bquote("Empirical failure intensity vs. cumulative test time:" ~.(SheetName))
)
```

\newpage
## Laplace Trend Test
`r if (x==1){paste0("The following figure shows the Laplace test statistic for reliability growth as a function of cumulative test time (FT). A decreasing trend indicates reliability growth, while an increasing trend indicates reliabilty deterioration. The Laplace test statistic on the y-axis corresponds to the critical values of a normal distribution. This means that if the trend falls below a specific level, then we cannot reject the null hypothesis that the failure data suggests reliability growth at a specified level of confidence. The six black dot-dash style lines correspond to the 90%, 95%, 99%, 99.9%, 99.9999%, and 99.999999% respectively. The red line is user-specified and has been set to ",percent_data_for_PSSE*100,"%. The level of confidence is a subjective choice made by the analyst. Reliability growth is desired because software reliability growth models assume curves that exhibit increasing time between failures. If reliability growth is not present than the model fitting step may fail or produce predictions that are inaccurate. Therefore, the Laplace test statistic provides an objective quantitative measure for the analyst to decide if predictions may or may not be accurate.")}`

```{r, echo=FALSE}
laplace_trend_test <- function(inter_failure) {
  n <- length(inter_failure)
  failure_time <- IF_to_FT(inter_failure)   
  laplace_trend <- c()
  laplace_trend[1] <- 0
  for(i in 2:n) {
    sumint <- 0
    for(j in 1:(i-1)) {
      sumint <- sumint + failure_time[j]  
    }
    laplace_trend[i] <-(((1/(i-1))*sumint) -(failure_time[i]/2))/(failure_time[i]*(1/(12*(i-1))^(0.5)))
  }
  trend_data <- data.frame(c(1:length(laplace_trend)),laplace_trend)
  names(trend_data) <- c("Index","Laplace_factor")
  return(trend_data)
}
LTT <- laplace_trend_test(IF)

# Two Tailed test
Confidence.lvl <- confidence_lvl # <-- comes from report-specifications.R
Significance <- qnorm(1 - Confidence.lvl)
Z.value.1 <- qnorm(0.1)
Z.value.2 <- qnorm(0.05)
Z.value.3 <- qnorm(0.01)
Z.value.4 <- qnorm(0.001)
Z.value.5 <- qnorm(0.0000001)
Z.value.6 <- qnorm(0.000000001)

#Display significance lines
localTrendPlot <- 
  ggplot(data = LTT, aes(Index, Laplace_factor)) + 
  geom_step() +
  theme_classic() +
  theme(plot.title = element_text(hjust = 0.5), plot.subtitle = element_text(hjust = 0.5)) +
  xlab("Failure Number") +
  scale_x_continuous(breaks = round(seq(0 , max(LTT$Index), by = 20), 1)) +
  ylab("Laplace Test Statistic") +
  ggtitle(bquote("Laplace Trend Test" ~.(SheetName))) + 
  labs(subtitle = bquote("Confidence =" ~.(Confidence.lvl))) +
  geom_hline(yintercept = Significance, color="Red") +
  geom_hline(yintercept = Z.value.1, color="black", linetype = "dotdash", alpha=0.8) +
  geom_hline(yintercept = Z.value.2, color="black", linetype = "dotdash", alpha=0.8) +
  geom_hline(yintercept = Z.value.3, color="black", linetype = "dotdash", alpha=0.8) +
  geom_hline(yintercept = Z.value.4, color="black", linetype = "dotdash", alpha=0.8) +
  geom_hline(yintercept = Z.value.5, color="black", linetype = "dotdash", alpha=0.8) +
  geom_hline(yintercept = Z.value.6, color="black", linetype = "dotdash", alpha=0.8)

localTrendPlot
```

\newpage
## Running arithmetic average

`r if (x==1){paste0("The running arithmetic average plots the average of the first k times between failure. An increasing trends indicates reliability growth, while a decreasing trend indicates reliability deterioration. This is intuitive because if the time between failures is increasing, then later failures will increase the average.")}`

```{r,echo=FALSE}
running_average_test <- function(inter_failure) {
  n <- length(inter_failure)
  runningAverage <- c()
  for(i in 1:n) {
    sum1 <-0
    for(j in 1:i) {
      sum1 <- sum1 + inter_failure[j]
    }
    runningAverage[i] <- (1/i)*sum1;
  }
  runningAverage <- data.frame(c(1:length(runningAverage)),runningAverage)
  names(runningAverage) <- c("Index","Running_Average")
  return(runningAverage)
}
plot(running_average_test(IF), type="s",xlab="Failure number", ylab="Running average of interfailure times",main=bquote("Running arithmetic average test:"~.(SheetName)))
```

\newpage
# Tab2: Set Up and Apply Models

```{r, echo=FALSE}
source('./SFRAT/models/DSS/DSS_BM_FT.R')
source('./SFRAT/models/DSS/Model_specifications.R')
source('./SFRAT/models/GO/GO_BM_FT.R')
source('./SFRAT/models/GO/Model_specifications.R')
source('./SFRAT/models/JM/JM_BM.R')
source('./SFRAT/models/JM/Model_specifications.R')
source('./SFRAT/models/GM/GM_BM.R')
source('./SFRAT/models/GM/Model_specifications.R')
source('./SFRAT/models/Wei/Wei_NM_FT.R')
source('./SFRAT/models/Wei/Model_specifications.R')
source('./SFRAT/models/ISS/ISS_NM_FT.R')
source('./SFRAT/models/ISS/Model_specifications.R')
source('./SFRAT/utility/prediction/Detailed_prediction.R')
source('./SFRAT/utility/metrics/GOF.R')

models <- c() #used for legend
modelListMVF <- list() #passed to plot
modelListMTTF <- list()
modelListFI <- list()
modelListRel <- list()
ftList <- list()
modelIterator <- 1;
DSS = c()
if('DSS' %in% models_to_apply && DSS_BM_MLE(FT)!="nonconvergence"&&!is.na(DSS_lnL(FT, DSS_BM_MLE(FT)))&& !is.na(DSS_BM_MLE(FT))) { 
  pred_failure_times <- get_prediction_t("DSS",DSS_BM_MLE(FT),num_failures_future_prediction,last(FT),length(FT))
  pred_failure_times[is.na(pred_failure_times)] <- 0
  for(i in 1:num_failures_future_prediction){
    pred_failure_times[i] <- pred_failure_times[i] + last(FT)
  }
  complete_ft <- c(FT, pred_failure_times)
  temp_d <- data.frame("FT"=complete_ft)

  modelListMVF[[modelIterator]] = DSS_MVF(DSS_BM_MLE(FT),temp_d)[,2]
  modelListMTTF[[modelIterator]] = DSS_MTTF(DSS_BM_MLE(FT),temp_d)[,2]
  modelListFI[[modelIterator]] = DSS_FI(DSS_BM_MLE(FT),temp_d)[,2]
  modelListRel[[modelIterator]] = DSS_R_growth(DSS_BM_MLE(FT),temp_d,mission_time)[,2]
  
  ftList[[modelIterator]]=complete_ft
  models <- append(models, "DSS")
  modelIterator <- modelIterator + 1
}
GO = c()
if('GO' %in% models_to_apply && GO_BM_MLE(FT)!="nonconvergence"&&!is.na(GO_lnL(FT, GO_BM_MLE(FT)))&&!is.na(GO_BM_MLE(FT))) {
  pred_failure_times <- get_prediction_t("GO",GO_BM_MLE(FT),num_failures_future_prediction,last(FT),length(FT))
  pred_failure_times[is.na(pred_failure_times)] <- 0
  for(i in 1:num_failures_future_prediction){
    pred_failure_times[i] <- pred_failure_times[i] + last(FT)
  }
  complete_ft <- c(FT, pred_failure_times)
  temp_d <- data.frame("FT"=complete_ft)

  modelListMVF[[modelIterator]] = GO_MVF(GO_BM_MLE(FT),temp_d)[,1]
  modelListMTTF[[modelIterator]] = GO_MTTF(GO_BM_MLE(FT),temp_d)[,2]
  modelListFI[[modelIterator]] = GO_FI(GO_BM_MLE(FT),temp_d)[,2]
  modelListRel[[modelIterator]] = GO_R_growth(GO_BM_MLE(FT),temp_d,mission_time)[,2]
  
  ftList[[modelIterator]]=complete_ft
  models <- append(models, "GO")
  modelIterator <- modelIterator + 1
}
JM = c()
if('JM' %in% models_to_apply && JM_BM_MLE(IF)!="nonconvergence"&&!is.na(JM_lnL(IF, JM_BM_MLE(IF)))&&!is.na(JM_BM_MLE(IF))) { 
  pred_failure_times <- get_prediction_t("JM",JM_BM_MLE(IF),num_failures_future_prediction,last(FT),length(FT))
  pred_failure_times[is.na(pred_failure_times)] <- 0
  for(i in 1:num_failures_future_prediction){
    pred_failure_times[i] <- pred_failure_times[i] + last(FT)
  }
  complete_ft <- c(FT, pred_failure_times)
  complete_if <- FT_to_IF(complete_ft)
  temp_d <- data.frame("IF"=complete_if, "FT"=complete_ft)
  
  modelListMVF[[modelIterator]] = JM_MVF(JM_BM_MLE(IF),temp_d)[,1]
  modelListMTTF[[modelIterator]] = JM_MTTF(JM_BM_MLE(IF),temp_d)[,2]
  modelListFI[[modelIterator]] = JM_FI(JM_BM_MLE(IF),temp_d)[,2]
  modelListRel[[modelIterator]] = JM_R_growth(JM_BM_MLE(FT),temp_d,mission_time)[,2]
  
  ftList[[modelIterator]]=complete_ft
  models <- append(models, "JM")
  modelIterator <- modelIterator + 1
}
GM = c()
if('GM' %in% models_to_apply && GM_BM_MLE(IF)!="nonconvergence"&&!is.na(GM_lnL(IF, GM_BM_MLE(IF)))&&!is.na(GM_BM_MLE(IF))) { 
  pred_failure_times <- get_prediction_t("GM",GM_BM_MLE(IF),num_failures_future_prediction,last(FT),length(FT))
  pred_failure_times[is.na(pred_failure_times)] <- 0
  for(i in 1:num_failures_future_prediction){
    pred_failure_times[i] <- pred_failure_times[i] + last(FT)
  }
  complete_ft <- c(FT, pred_failure_times)
  complete_if <- FT_to_IF(complete_ft)
  temp_d <- data.frame("IF"=complete_if, "FT"=complete_ft)
  
  modelListMVF[[modelIterator]] = GM_MVF(GM_BM_MLE(IF),temp_d)[,1]
  modelListMTTF[[modelIterator]] = GM_MTTF(GM_BM_MLE(IF),temp_d)[,2]
  modelListFI[[modelIterator]] = GM_FI(GM_BM_MLE(IF),temp_d)[,2]
  modelListRel[[modelIterator]] = GM_R_growth(GM_BM_MLE(FT),temp_d,mission_time)[,2]
  
  ftList[[modelIterator]]=complete_ft
  models <- append(models, "GM")
  modelIterator <- modelIterator + 1
}
Wei = c()
if('Wei' %in% models_to_apply && Wei_NM_MLE(FT)!="nonconvergence"&&!is.na(Wei_lnL(FT, Wei_NM_MLE(FT)))&&!is.na(Wei_NM_MLE(FT))) { 
  pred_failure_times <- get_prediction_t("Wei",Wei_NM_MLE(FT),num_failures_future_prediction,last(FT),length(FT))
  pred_failure_times[is.na(pred_failure_times)] <- 0
  for(i in 1:num_failures_future_prediction){
    pred_failure_times[i] <- pred_failure_times[i] + last(FT)
  }
  complete_ft <- c(FT, pred_failure_times)
  temp_d <- data.frame("FT"=complete_ft)
  
  modelListMVF[[modelIterator]] = Wei_MVF(Wei_NM_MLE(FT),temp_d)[,1]
  modelListMTTF[[modelIterator]] = Wei_MTTF(Wei_NM_MLE(FT),temp_d)[,2]
  modelListFI[[modelIterator]] = Wei_FI(Wei_NM_MLE(FT),temp_d)[,2]
  modelListRel[[modelIterator]] = Wei_R_growth(Wei_NM_MLE(FT),temp_d,mission_time)[,2]
  
  ftList[[modelIterator]]=complete_ft
  models <- append(models, "Wei") 
  modelIterator <- modelIterator + 1
}

ISS = c()
if('ISS' %in% models_to_apply && ISS_NM_MLE(FT)!="nonconvergence"&&!is.na(ISS_lnL(FT, ISS_NM_MLE(FT)))&&!is.na(ISS_NM_MLE(FT))) { 
  pred_failure_times <- get_prediction_t("ISS",ISS_NM_MLE(FT),num_failures_future_prediction,last(FT),length(FT))
  pred_failure_times[is.na(pred_failure_times)] <- 0
  for(i in 1:num_failures_future_prediction){
    pred_failure_times[i] <- pred_failure_times[i] + last(FT)
  }
  complete_ft <- c(FT, pred_failure_times)
  temp_d <- data.frame("FT"=complete_ft)
  
  modelListMVF[[modelIterator]] = ISS_MVF(ISS_NM_MLE(FT),temp_d)[,1]
  modelListMTTF[[modelIterator]] = ISS_MTTF(ISS_NM_MLE(FT),temp_d)[,2]
  modelListFI[[modelIterator]] = ISS_FI(ISS_NM_MLE(FT),temp_d)[,2]
  modelListRel[[modelIterator]] = ISS_R_growth(ISS_NM_MLE(FT),temp_d,mission_time)[,2]
  
  ftList[[modelIterator]]=complete_ft
  models <- append(models, "ISS") 
  modelIterator <- modelIterator + 1
}

hxyMVF <- data.frame(modelListMVF)
hxyMTTF <- data.frame(modelListMTTF)
hxyMTTF[is.na(hxyMTTF)] <- 0
hxyFI <- data.frame(modelListFI)
hxyFI[is.na(hxyFI)] <- 0
hxyRel <- data.frame(modelListRel)
hxyRel[is.na(hxyRel)] <- 0
hft <- data.frame(ftList)
ntrees <- length(models)
```

## Cumulative failures
`r if (x==1){paste0("The following figure shows the fit of",if('DSS' %in% models){paste(" delayed s-shaped, ")}, if('GM' %in% models){paste(" geometric, ")},if('Wei' %in% models){paste(" Weibull, ")},if('GO' %in% models){paste(" Goel-Okumoto, ")},if('JM' %in% models){paste(" Jelinski-Moranda,  ")},if('ISS' %in% models){paste(" inflexion s-shaped")}, "models to the cumulative number of failures detected in the ",SheetName," data.")}`

```{r, echo=FALSE}
plot(hft[,which(hft==max(hft))/dim(hft)[1]],1:dim(hft)[1], type="n",xlab="Cumulative test time", ylab="Cumulative failures")
linetype <- c("p","l","o","b","c","s","S","h")  
lines(FT, FN, type="s")
for (i in 1:ntrees) {
  lines(hft[,i], hxyMVF[,i], type="l", col=colors[i])
}
abline(v = tail(FT, n=1))
title(bquote("Cumulative failures vs. cumulative test time:" ~.(SheetName)))
legend("topleft", legend=models, lty=c(1,1), cex=0.8, col=colors)
```

\newpage
## Times between failures

`r if (x==1){paste0("The following figure shows the fit of",if('DSS' %in% models){paste(" delayed s-shaped, ")}, if('GM' %in% models){paste(" geometric, ")},if('Wei' %in% models){paste(" Weibull, ")},if('GO' %in% models){paste(" Goel-Okumoto, ")},if('JM' %in% models){paste(" Jelinski-Moranda,  ")},if('ISS' %in% models){paste(" inflexion s-shaped")}, "models to the times between failures detected in the ",SheetName," data.")}`

```{r,echo=FALSE}
#Plot 2
hif <- data.frame(FT_to_IF(complete_ft)) 
plot(hft[,which(hft==max(hft))/dim(hft)[1]],FT_to_IF(hft[,which(hft==max(hft))/dim(hft)[1]]), type="n",xlab="Cumulative test time", ylab="Times between failures")
lines(FT, IF, type="s")
for (i in 1:ntrees) {
  lines(hft[,i], hxyMTTF[,i], type="l", col=colors[i])
}
abline(v = tail(FT, n=1))
title(bquote("Times between failures vs. cumulative test time:" ~.(SheetName)))
legend("topleft", legend=models, lty=c(1,1), cex=0.8, col=colors)
```

\newpage
## Failure intensity

`r if (x==1){paste0("The following figure shows the fit of",if('DSS' %in% models){paste(" delayed s-shaped, ")}, if('GM' %in% models){paste(" geometric, ")},if('Wei' %in% models){paste(" Weibull, ")},if('GO' %in% models){paste(" Goel-Okumoto, ")},if('JM' %in% models){paste(" Jelinski-Moranda,  ")},,if('ISS' %in% models){paste(" inflexion s-shaped")}, "models to the failure intensity of the ",SheetName," data.")}`

```{r, echo=FALSE}
tempFI <- 1/(FT_to_IF(hft[,which(hft==max(hft))/dim(hft)[1]]))
tempFI[is.infinite(tempFI)] <- 0
upLimTemp <- list()
for(i in 1:ntrees){
  upLimTemp <- max(hxyFI[,i])
}

upLim <- max(upLimTemp)
plot(
  hft[,which(hft==max(hft))/dim(hft)[1]], tempFI, type="n",ylim=c(min(tempFI),min(upLim,max(tempFI))),
  xlab="Cumulative test time", ylab="Failure intensity"
)
lines(FT, 1/IF, type="s")
for (i in 1:ntrees) {
  lines(hft[,i], hxyFI[,i], type="l", col=colors[i])
}
abline(v = tail(FT, n=1))
title(bquote("Failure intensity vs. Cumulative test time" ~.(SheetName)))
legend("topright", legend=models, lty=c(1,1), cex=0.8, col=colors)
```

\newpage
## Reliability growth


`r if (x==1){paste0("The following figure shows the reliability growth curve of the  fit of",if('DSS' %in% models){paste(" delayed s-shaped, ")}, if('GM' %in% models){paste(" geometric, ")},if('Wei' %in% models){paste(" Weibull, ")},if('GO' %in% models){paste(" Goel-Okumoto, ")},if('JM' %in% models){paste(" Jelinski-Moranda,")},if('ISS' %in% models){paste(" inflexion s-shaped")}, "models to the ",SheetName," data. The data itself does not display. This plot indicates a models prediction that the software will be reliable (exhibit zero failures) for a duration of ",mission_time," time units as a function of cumulative test time (FT). Selecting a model upon which to base a reliability assessment is a subjective choice made by the analyst. Statistical measures of goodness of fit, reported on page 12 of this report can be used to decide this decision making process. If the Laplace test statistic does not exhibit reliability growth, than a conservative approach is to document this as the reason why no reliability estimate is provided at the time of preparing a report.")}`

```{r,echo=FALSE}
plot(hft[,which(hft==max(hft))/dim(hft)[1]], seq(from=0, to=max(hxyRel),by=max(hxyRel)/(length(complete_ft)-1)),type="n", xlab="Cumulative test time", ylab="Reliability growth")
for (i in 1:ntrees) { 
  lines(hft[,i], hxyRel[,i], type="l", col=colors[i]) ##Limit the upperbound to maximum value
} 
abline(v = tail(FT, n=1))
title(bquote("Reliability growth vs. cumulative test time" ~.(SheetName)))
legend("topleft", legend=models, lty=c(1,1), cex=0.8, col=colors)
```

\newpage
# Tab3: Query Model Results
`r if (x==1){paste0("The following table shows inferences enabled by the models, including the time to achieve a reliability of ",desired_reliability*100,"% (probability of zero failures for ",reliability_interval_length," time units), expected number of failures in the next ",additional_time_software_will_run," time units, and expected time to observe an additional ",num_failures_to_predict," failures computed for the fit of",if('DSS' %in% models){paste(" delayed s-shaped (DSS), ")}, if('GM' %in% models){paste(" geometric (GM), ")},if('Wei' %in% models){paste(" Weibull (Wei), ")},if('GO' %in% models){paste(" Goel-Okumoto (GO), ")},if('JM' %in% models){paste(" Jelinski-Moranda (JM), ")},if('ISS' %in% models){paste(" inflexion s-shaped (ISS) ")}, "models to the models.")}`

```{r,echo=FALSE}
predTable <- matrix(NA, nrow = length(models) , ncol = 3) #Setting Up The Table
rownames(predTable) <- models
colnames(predTable) <- c("Time to achieve specified reliability","Expected number of failures","Expected time to N failure")
options(digits=4)   ##Calulating Future Failures #Time to achieve specified reliability

firstCol <- c()     #First column of prediction table
secondCol <- c() #Expected number of failures
thirdCol <- c()#Expected time to N failure
if('DSS' %in% models) { 
  firstCol=append(firstCol,get_reliability_t("DSS",DSS_BM_MLE(FT),desired_reliability,reliability_interval_length,last(FT),num_failures_to_predict))
  secondCol=append(secondCol, get_prediction_k("DSS",DSS_BM_MLE(FT),additional_time_software_will_run,last(FT),length(FT)))
  thirdCol=append(thirdCol,last(get_prediction_t("DSS",DSS_BM_MLE(FT),num_failures_to_predict,last(FT),length(FT))))
  }
if('GO' %in% models) {
  firstCol=append(firstCol,get_reliability_t("GO",GO_BM_MLE(FT),desired_reliability,reliability_interval_length,last(FT),num_failures_to_predict))
  secondCol=append(secondCol,get_prediction_k("GO",GO_BM_MLE(FT),additional_time_software_will_run,last(FT),length(FT)))
  thirdCol=append(thirdCol,last(get_prediction_t("GO",GO_BM_MLE(FT),num_failures_to_predict,last(FT),length(FT))))
}
if('JM' %in% models) { 
  firstCol=append(firstCol,get_reliability_t("JM",JM_BM_MLE(IF),desired_reliability,reliability_interval_length,last(IF),num_failures_to_predict))
  secondCol=append(secondCol,get_prediction_k("JM",JM_BM_MLE(IF),additional_time_software_will_run,last(FT),length(FT)))
  thirdCol=append(thirdCol,last(get_prediction_t("JM",JM_BM_MLE(IF),num_failures_to_predict,last(FT),length(FT))))
}
if('GM' %in% models) {
  firstCol=append(firstCol,get_reliability_t("GM",GM_BM_MLE(IF),desired_reliability,reliability_interval_length,last(FT),num_failures_to_predict))##The right hand side is creating a null model. TODO: Refer back to the original file to solve this precision issue
  secondCol=append(secondCol,get_prediction_k("GM",GM_BM_MLE(IF),additional_time_software_will_run,last(FT),length(FT)))
  thirdCol=append(thirdCol,last(get_prediction_t("GM",GM_BM_MLE(IF),num_failures_to_predict,last(FT),length(FT))))
}
if('Wei' %in% models) { 
  firstCol=append(firstCol,get_reliability_t("Wei",Wei_NM_MLE(FT),desired_reliability,reliability_interval_length,last(FT),num_failures_to_predict))
  secondCol=append(secondCol, get_prediction_k("Wei",Wei_NM_MLE(FT),additional_time_software_will_run,last(FT),length(FT)))
  thirdCol=append(thirdCol,last(get_prediction_t("Wei",Wei_NM_MLE(FT),num_failures_to_predict,last(FT),length(FT))))
}

if('ISS' %in% models) { 
  firstCol=append(firstCol,get_reliability_t("ISS",ISS_NM_MLE(FT),desired_reliability,reliability_interval_length,last(FT),num_failures_to_predict))
  secondCol=append(secondCol, get_prediction_k("ISS",ISS_NM_MLE(FT),additional_time_software_will_run,last(FT),length(FT)))
  thirdCol=append(thirdCol,last(get_prediction_t("ISS",ISS_NM_MLE(FT),num_failures_to_predict,last(FT),length(FT))))
}



predTable[,1] <- firstCol
predTable[,2] <- secondCol
predTable[,3] <- thirdCol
kable(predTable,align='r',digits=4)
```

\newpage
# Tab4: Evaluate Models

`r if (x==1){paste0("The following table shows the measures of goodness of fit computed for the ", if('DSS' %in% models){paste(" delayed s-shaped, ")}, if('GM' %in% models){paste(" geometric, ")},if('Wei' %in% models){paste(" Weibull, ")},if('GO' %in% models){paste(" Goel-Okumoto, ")},if('JM' %in% models){paste(" Jelinski-Moranda, ")}, if('ISS' %in% models){paste(" inflexion s-shaped ")}, " The Akaike Information Criterion (AIC) is an information theoretic measure. Lower values are preferred. The GM model achieved the lowest AIC value on the ",SheetName," data. A difference of 2.0 or more in the AIC values of two models indicates the model with the lower AIC score is preferred with statistical significance. The Predictive Sum of Squares Error (PSSE) used ",percent_data_for_PSSE*100,"% of the ",SheetName," data to fit the models and computed the sum of the squares between the differences of the remaining 10% of the data not used to fit the models. Lower values are preferred. The GM model achieved the lowest PSSE value on the ",SheetName," data. The measures of goodness of fit can help select a model, but the choice is ultimately a subjective choice made by the analyst.")}`

```{r,echo=FALSE}
GOFTable <- matrix(NA, nrow = length(models), ncol = 2)#Setup Table
rownames(GOFTable) <- models
colnames(GOFTable) <- c(bquote("Akaike Information Criterion (AIC)"), bquote("Predictive sum of squares error (PSSE)" ~.(percent_data_for_PSSE)))

GOFFirstCol <- c()#Calculate AIC
GOFSecondCol <- c() #Calculate PSSE
if('DSS' %in% models) { 
  GOFFirstCol=append(GOFFirstCol,aic(2, DSS_lnL(FT, DSS_BM_MLE(FT))))
  GOFSecondCol=append(GOFSecondCol,psse("DSS", FT, DSS_BM_MLE(FT), percent_data_for_PSSE))
}
if('GO' %in% models) {
  GOFFirstCol=append(GOFFirstCol,aic(2, GO_lnL(FT, GO_BM_MLE(FT))))
  GOFSecondCol=append(GOFSecondCol,psse("GO", FT, GO_BM_MLE(FT), percent_data_for_PSSE))
}
if('JM' %in% models) { 
  GOFFirstCol=append(GOFFirstCol,aic(2, JM_lnL(IF, JM_BM_MLE(IF))))
  GOFSecondCol=append(GOFSecondCol, psse("JM", FT, JM_BM_MLE(IF), percent_data_for_PSSE))
}
if('GM' %in% models) {
  GOFFirstCol=append(GOFFirstCol,aic(2, GM_lnL(IF, GM_BM_MLE(IF))))
  GOFSecondCol=append(GOFSecondCol,psse("GM", FT, GM_BM_MLE(IF), percent_data_for_PSSE))
}
if('Wei' %in% models) { 
  GOFFirstCol=append(GOFFirstCol,aic(3, Wei_lnL(FT, Wei_NM_MLE(FT))))
  GOFSecondCol=append(GOFSecondCol, psse("Wei", FT, Wei_NM_MLE(FT), percent_data_for_PSSE))
}
if('ISS' %in% models) { 
  GOFFirstCol=append(GOFFirstCol,aic(3, ISS_lnL(FT, ISS_NM_MLE(FT))))
  GOFSecondCol=append(GOFSecondCol, psse("ISS", FT, ISS_NM_MLE(FT), percent_data_for_PSSE))
}

GOFFirstCol <-  round(GOFFirstCol,digits=2)
fircol <- c()
for(i in 1:length(models)){
  if(GOFFirstCol[i]==min(GOFFirstCol)){
    fircol <- append(fircol,paste(c("*",GOFFirstCol[i]),collapse=""))
  }else{
    fircol <- append(fircol,GOFFirstCol[i])
  }
  i=i+1;
}
GOFTable[,1] <- fircol

GOFSecondCol <-  round(GOFSecondCol,digits=2)
seccol <- c()
for(i in 1:length(models)){
  if(GOFSecondCol[i]==min(GOFSecondCol)){
    seccol <- append(seccol,paste(c("*",GOFSecondCol[i]),collapse=""))
  }else{
    seccol <- append(seccol,GOFSecondCol[i])
  }
  i=i+1;
}
GOFTable[,2] <- seccol

kable(GOFTable,align='r')
```