diff --git a/R-package/R/inner_procedure.r b/R-package/R/inner_procedure.r
index bb804538e7b50a73f32c172635ef581e8b8b10df..bc702c306dbff3808dc7f5378fb26e8d8d352df7 100644
--- a/R-package/R/inner_procedure.r
+++ b/R-package/R/inner_procedure.r
@@ -20,11 +20,12 @@
 #' @param max_num_cores maximum number of cores for parallelism
 #' @param progress_bar if true use a normal progress bar, otherwise a shiny's progress bar
 #'
-#' @return a list with 4 elements:
+#' @return a list with 5 elements:
 #' 1) returns 0 succeeded, 1 increase init.samples, or 2 increase init.samples and nb_models
 #' 2) An index
-#' 3) Accuracy vector
-#' 4) A vector with a Habitat objects, each consisting of 7 slots: \cr
+#' 3) num_models number of models selected
+#' 4) Accuracy vector
+#' 5) A vector with a Habitat objects, each consisting of 7 slots: \cr
 #' run1@models - list of selected classifiers (only if save_runs is TRUE) \cr
 #' run1@ref_samples - list of SpatialPointsDataFrames with same length as run1@models holding reference labels [1,2] for each selected model (only if save_runs is TRUE) \cr
 #' run1@switch - vector of length run1@models indicating if target class equals 2, if not NA the labels need to be switched (only if save_runs is TRUE) \cr
@@ -168,7 +169,13 @@ sample_nb <- function(raster,
         length(which_models_null[which_models_null == FALSE]) == length(models_list)) {
       remove(points_list)
       remove(models_list)
-      out <- list(returns = 1, index = NULL, acc = NULL, obj = NULL)
+      out <- list(
+        returns = 1,
+        index = NULL,
+        num_models = 0,
+        acc = NULL,
+        obj = NULL,
+      )
       return(out)
     }
     if (length(which_models_null[which_models_null == FALSE]) > 0) {
@@ -228,11 +235,18 @@ sample_nb <- function(raster,
     index <- which.max(dif[2,])
     ch <- as.numeric(na.omit(channel[, index]))
     if (length(ch) == 0) {
-      out <- list(returns = 2, index = NULL, acc = NULL, obj = NULL)
+      out <- list(
+        returns = 2,
+        index = NULL,
+        num_models = 0,
+        acc = NULL,
+        obj = NULL
+      )
       return(out)
     }
     acc <- (round(m[l] ^ 2, 2) / 0.25)
 
+    cat("\n")
     print(paste("class=", index, "  difference=", (round(m[l] ^ 2, 2) / 0.25),
                 sep = ""))
     l <- l + 1
@@ -248,6 +262,7 @@ sample_nb <- function(raster,
   }
   models <- models[ch]
   print(paste("n_models =", length(models)))
+  flush(stdout())
   switch <- switch[ch, index]
   points <- points_list[ch]
   remove(points_list)
@@ -275,7 +290,7 @@ sample_nb <- function(raster,
           raster::reclassify(result1[[j]], rbind(c(0.5, 1.5, 2), c(1.6, 2.5, 1)))
       }
     }
-    print(j)
+    #print(j)
     j <- j + 1
   }
   ###
@@ -306,6 +321,7 @@ sample_nb <- function(raster,
       seeds = 0
     )
   }
+  num_models <- length(models)
   remove(models)
   remove(points)
   remove(switch)
@@ -315,6 +331,12 @@ sample_nb <- function(raster,
   remove(seed2)
   gc(full = TRUE)
 
-  out <- list(returns = 0, index = index, acc = acc, obj = obj)
+  out <- list(
+    returns = 0,
+    index = index,
+    num_models = num_models,
+    acc = acc,
+    obj = obj
+  )
   return(out)
 }
diff --git a/R-package/R/outer_procedure.r b/R-package/R/outer_procedure.r
index 083d89d2fea622ee51ac07dc6967d194fc2de3e8..1b9d87c66c4846c63fbf9fbdd227d1314d933bea 100644
--- a/R-package/R/outer_procedure.r
+++ b/R-package/R/outer_procedure.r
@@ -125,6 +125,7 @@ multi_Class_Sampling <- function(in.raster,
     }
   }
 
+  print(paste(paste("Habitat", 0), "Starting"))
   for (i in step:r) {
     if (i == r) {
       last = T
@@ -164,6 +165,7 @@ multi_Class_Sampling <- function(in.raster,
             )
             returns <- maFo_rf$returns
             index <- maFo_rf$index
+            num_models <- maFo_rf$num_models
             acc <- maFo_rf$acc
             maFo_rf <- maFo_rf$obj
             if (returns == 1) {
@@ -263,6 +265,7 @@ multi_Class_Sampling <- function(in.raster,
         )
         returns <- maFo_rf$returns
         index <- maFo_rf$index
+        num_models <- maFo_rf$num_models
         acc <- maFo_rf$acc
         maFo_rf <- maFo_rf$obj
         if (returns == 1) {
@@ -318,6 +321,7 @@ multi_Class_Sampling <- function(in.raster,
       r = RGB[1],
       g = RGB[2],
       b = RGB[3],
+      num_models = num_models,
       acc = acc,
       color = color,
       outPath = outPath,
@@ -369,6 +373,7 @@ multi_Class_Sampling <- function(in.raster,
         )
         returns <- maFo_rf$returns
         index <- maFo_rf$index
+        num_models <- maFo_rf$num_models
         acc <- maFo_rf$acc
         maFo_rf <- maFo_rf$obj
         if (returns == 1) {
@@ -423,6 +428,7 @@ multi_Class_Sampling <- function(in.raster,
         r = RGB[1],
         g = RGB[2],
         b = RGB[3],
+        num_models = num_models,
         acc = acc,
         color = color,
         outPath = outPath,
@@ -470,6 +476,8 @@ multi_Class_Sampling <- function(in.raster,
 
     thres <- as.numeric(decision)
     dummy <- maFo_rf@layer[[1]]
+    max_prob <- raster::cellStats(dummy, "max")
+    thres <- (max_prob/100) * thres
     dummy[dummy < thres] <- 1
     dummy[dummy >= thres] <- NA
     in.raster <- in.raster * dummy
@@ -505,6 +513,9 @@ multi_Class_Sampling <- function(in.raster,
       break()
     }
 
+    num_habitat <- i + 1
+    print(paste("Habitat", num_habitat, "Starting", sep = " "))
+    flush(stdout())
     decision2 <- readline(paste("Adjust init.samples/nb.models (actual ",
                                 init.samples,
                                 "/",
diff --git a/R-package/R/plot_interactive.r b/R-package/R/plot_interactive.r
index 8dac47a5b17d72e217a086f475e01cf2cd31d8f5..9dbf62fc00b71b7afcd1dbcb9ff59ef4a9713fa4 100644
--- a/R-package/R/plot_interactive.r
+++ b/R-package/R/plot_interactive.r
@@ -9,13 +9,25 @@
 #' @param g green channel (integer)
 #' @param b blue channel (integer)
 #' @param acc predictive accuracy (integer)
+#' @param num_models number of selected models
 #' @param color color pallet
 #' @param outPath file path for '.html export (character)
 #' @param plot_on_browser plot on the browser or inline in a notebook (default TRUE)
 #' @return a leafLet map
 #'
 #' @export
-iplot <- function(x, y, HaTy, r, g, b, acc, color, outPath, plot_on_browser = TRUE) {
+iplot <- function(x,
+                  y,
+                  HaTy,
+                  r,
+                  g,
+                  b,
+                  acc,
+                  num_models,
+                  color,
+                  outPath,
+                  plot_on_browser = TRUE) {
+  
   #x=layerInfo, y=RGB Image
   ##############################################################################
   if (exists("color") == F) {
@@ -113,7 +125,9 @@ iplot <- function(x, y, HaTy, r, g, b, acc, color, outPath, plot_on_browser = TR
   ######
   rr <- x
   raster::values(rr) <- 1:raster::ncell(rr)
-
+  min <- raster::cellStats(x, "min")
+  max <- raster::cellStats(x, "max")
+  x[x<min] <- NA
   ##############################################################################
   ##[2] Create Leaflet Html output for Webbrowser
   mv <- leaflet::leaflet(options = leaflet::leafletOptions(zoomControl = FALSE))
@@ -154,8 +168,8 @@ iplot <- function(x, y, HaTy, r, g, b, acc, color, outPath, plot_on_browser = TR
     map = mv,
     "bottomright",
     pal = pal,
-    values = raster::cellStats(x, "range"),
-    title = "Habitat Type Probability",
+    values =  c((min*100)/max, 100),
+    title = paste("Habitat Type Probability<br>(for ", num_models, " models)", sep = ""),
     opacity = 1
   )
   mv <- leaflet::addLayersControl(map = mv,
diff --git a/R-package/man/iplot.Rd b/R-package/man/iplot.Rd
index 7f48ffe8bda0725e2bffd55e056344e6a12f884c..d0bdd587e9ae52efbb3240557bc23d56530f3d37 100644
--- a/R-package/man/iplot.Rd
+++ b/R-package/man/iplot.Rd
@@ -4,7 +4,19 @@
 \alias{iplot}
 \title{Plot Habitat Types}
 \usage{
-iplot(x, y, HaTy, r, g, b, acc, color, outPath, plot_on_browser = TRUE)
+iplot(
+  x,
+  y,
+  HaTy,
+  r,
+  g,
+  b,
+  acc,
+  num_models,
+  color,
+  outPath,
+  plot_on_browser = TRUE
+)
 }
 \arguments{
 \item{x}{probability image (*rasterObject)}
@@ -21,6 +33,8 @@ iplot(x, y, HaTy, r, g, b, acc, color, outPath, plot_on_browser = TRUE)
 
 \item{acc}{predictive accuracy (integer)}
 
+\item{num_models}{number of selected models}
+
 \item{color}{color pallet}
 
 \item{outPath}{file path for '.html export (character)}
diff --git a/R-package/man/sample_nb.Rd b/R-package/man/sample_nb.Rd
index ce0baa5b8097d2b55f7fce8342f6f07732654633..171f2e4bd6a3273b2442e8c9b84b6d55a134c451 100644
--- a/R-package/man/sample_nb.Rd
+++ b/R-package/man/sample_nb.Rd
@@ -60,10 +60,11 @@ sample_nb(
 \item{nb_models}{number of models (independent classifiers) to collect}
 }
 \value{
-a list with 4 elements:
+a list with 5 elements:
 \enumerate{
 \item returns 0 succeeded, 1 increase init.samples, or 2 increase init.samples and nb_models
 \item An index
+\item num_models number of models selected
 \item Accuracy vector
 \item A vector with a Habitat objects, each consisting of 7 slots: \cr
 run1@models - list of selected classifiers (only if save_runs is TRUE) \cr