add option extract output that was not used for calibration, by assigning all observed data of that variable as NA

Author: nguyenta

R/displayOutput.R

@@ -214,24 +214,43 @@ printVariableNameObservedFiles <- function(outputExtraction){
 plotSimulated <- function(inputDataFrame){
   inputDataFrame$date <- as.Date(inputDataFrame$date, "%Y-%m-%d")
 
-  myplot <- ggplot(inputDataFrame) + 
-    #95 PPU
-    geom_ribbon(aes(x = date, ymin = lower, ymax = upper, color = "95PPU"), 
-                size = 0.0, fill = "red", alpha = 0.3) +
-    # Best simulation
-    geom_line(aes(x = date, y = best, color = "Best simulation"), alpha = 0.6) +
-    # Observed data
-    geom_line(aes(x = date, y = median, color = "Median"), alpha = 0.6) +
-    # Observed data
-    geom_line(aes(x = date, y = observed, color = "Observed"), alpha = 0.6) +
-    scale_colour_manual(name= '', values=c("95PPU" = "red", 
-                                           "Best simulation" = "red",
-                                           "Median" = "green",
-                                           "Observed" = "darkblue")) +
-    # Axis name
-    labs(x ="  ", y = " ") +
-    scale_x_date(date_labels = "%m-%Y") +
-    theme_bw()
+  if(length(which(is.na(inputDataFrame$observed))) == nrow(inputDataFrame)){
+    myplot <- ggplot(inputDataFrame) + 
+      #95 PPU
+      geom_ribbon(aes(x = date, ymin = lower, ymax = upper, color = "95PPU"), 
+                  size = 0.0, fill = "red", alpha = 0.3) +
+      # Best simulation
+      geom_line(aes(x = date, y = best, color = "Best simulation"), alpha = 0.6) +
+      # Observed data
+      geom_line(aes(x = date, y = median, color = "Median"), alpha = 0.6) +
+      scale_colour_manual(name= '', values=c("95PPU" = "red", 
+                                             "Best simulation" = "red",
+                                             "Median" = "green")) +
+      # Axis name
+      labs(x ="  ", y = " ") +
+      scale_x_date(date_labels = "%m-%Y") +
+      theme_bw()    
+  } else {
+    myplot <- ggplot(inputDataFrame) + 
+      #95 PPU
+      geom_ribbon(aes(x = date, ymin = lower, ymax = upper, color = "95PPU"), 
+                  size = 0.0, fill = "red", alpha = 0.3) +
+      # Best simulation
+      geom_line(aes(x = date, y = best, color = "Best simulation"), alpha = 0.6) +
+      # Observed data
+      geom_line(aes(x = date, y = median, color = "Median"), alpha = 0.6) +
+      # Observed data
+      geom_line(aes(x = date, y = observed, color = "Observed"), alpha = 0.6) +
+      scale_colour_manual(name= '', values=c("95PPU" = "red", 
+                                             "Best simulation" = "red",
+                                             "Median" = "green",
+                                             "Observed" = "darkblue")) +
+      # Axis name
+      labs(x ="  ", y = " ") +
+      scale_x_date(date_labels = "%m-%Y") +
+      theme_bw()    
+  }
+
   return(myplot)
 }
 

R/postProcessing.R

@@ -32,8 +32,8 @@ behaSimulation <- function(objValue, simData, parameterValue, behThreshold,
   ppuSimData <- matrix(rep(NA, 4*nrow), ncol = 4)
 
   if (samplingApproach == 'Cali_(Generalized_Likelihood_Uncertainty_Estimation)'){
-    # If GLUE then use different approach
-    # Normalized the likelikhood
+
+    # If GLUE then use different approach, normalized the likelikhood
     normalizedLik <- which(objValue >= behThreshold)
     normalizedLik <- normalizedLik/sum(normalizedLik)
     for (i in 1:nrow){
@@ -104,30 +104,38 @@ behaSimulation <- function(objValue, simData, parameterValue, behThreshold,
 prFactor <- function(obs, low, up){
 
   naIndex <- which(is.na(obs))
+  
   if (length(naIndex) > 0){
     obs <- obs[-c(naIndex)]
     low <- low[-c(naIndex)]
     up <- up[-c(naIndex)]    
   }
 
-  pfactor <- 0
-  count <- 0
-  count_all <- 0
-  rfactor <- 0
-  
-  for (i in 1:length(obs)){
-    count_all <- count_all + 1
-    rfactor <- rfactor + up[i] - low[i]
-    if(obs[i] >= low[i]){
-      if(obs[i] <= up[i]){
-        count <- count + 1
+  #If all no single value in the observed data
+  if (length(obs) == 0){
+    pfactor <- NA
+    rfactor <- NA
+    
+  } else {
+    pfactor <- 0
+    count <- 0
+    count_all <- 0
+    rfactor <- 0
+    
+    for (i in 1:length(obs)){
+      count_all <- count_all + 1
+      rfactor <- rfactor + up[i] - low[i]
+      if(obs[i] >= low[i]){
+        if(obs[i] <= up[i]){
+          count <- count + 1
+        }
       }
     }
+    
+    pfactor <- count/count_all
+    rfactor <- rfactor/(count_all * sd(obs))    
   }
 
-  pfactor = count/count_all
-  rfactor <- rfactor/(count_all * sd(obs))
-  
   return(c(pfactor, rfactor))
 }
 
@@ -233,7 +241,11 @@ calObjFunction <- function(parameterValue, ncores,
           simData[[j]][[counter[j]]] <- tempSimData[(sIndex + 1):eIndex, 1]
           output$simData[[j]][[counter[j]]] <- simData[[j]][[counter[j]]]
 
-          if (index != 'userObjFunction'){
+          # check if observed variables has no data, all are missing
+          missingValue <- which(is.na(observedData[[j]][,2]))
+          
+          if ((index != 'userObjFunction') & 
+              (length(missingValue) != length(observedData[[j]][,2]))){
 
             output$perCriteria[[j]][[counter[j]]] <- perCriteria(observedData[[j]][,2],
                                                                  simData[[j]][[counter[j]]])

data/parameter_range.csv

@@ -1228,11 +1228,14 @@ print(sensCaliObject)[]")
 
             } else {
 
-              if (is.numeric(temp$V2)){
+              if (is.numeric(temp$V2) |
+                  length(which(is.na(temp$V2))) == length(temp$V2)){
+                
                 globalVariable$observedData[[i]] <<- temp
                 colnames(globalVariable$observedData[[i]]) <<- c("Date", "Value")
 
               } else {
+                
                 checkGetObservedDataFileMessage <- paste("Error: Data in the second column are not numeric values: ",
                                                          globalVariable$observedDataFile[i],
                                                          sep = "")
@@ -1596,7 +1599,7 @@ print(sensCaliObject)[]")
       # Give column names for the 95PPU table
       colnames(tempVar) <- c("date", "lower", "median", "upper", "best", "observed")
       globalVariable$PlotVariableNumber <<- plotSimulated(tempVar)
-
+      
       # Plot the 95PPU
       output$PlotVariableNumber <- renderPlotly(
         ggplotly(globalVariable$PlotVariableNumber +