Diff of /R/ich_predict.R [000000] .. [242173]

Switch to side-by-side view

--- a
+++ b/R/ich_predict.R
@@ -0,0 +1,203 @@
+#' @title Predict ICH Images
+#' @description This function will take the data.frame of predictors and
+#' predict the ICH voxels from the model chosen.
+#'
+#' @param df \code{\link{data.frame}} of predictors.  If \code{multiplier}
+#' column does not exist, then \code{\link{ich_candidate_voxels}} will
+#' be called
+#' @param nim object of class \code{\link{nifti}}, from
+#' \code{\link{make_predictors}}
+#' @param model model to use for prediction,
+#' either the random forest (rf) or logistic
+#' @param verbose Print diagnostic output
+#' @param native Should native-space predictions be given?
+#' @param native_img object of class \code{\link{nifti}}, which
+#' is the dimensions of the native image
+#' @param transformlist Transforms list for the transformations back to native space.
+#' NOTE: these will be inverted.
+#' @param interpolator Interpolator for the transformation back to native space
+#' @param native_thresh Threshold for re-thresholding binary mask after
+#' interpolation
+#' @param shiny Should shiny progress be called?
+#' @param model_list list of model objects, used mainly for retraining
+#' but only expert use.
+#' @param smoothed_cutoffs A list with an element
+#' \code{mod.dice.coef}, only expert use.
+#' @param outfile filename for output file.
+#' We write the smoothed, thresholded image.  If \code{native = TRUE},
+#' then the file will be native space, otherwise in registered
+#' space
+#' @param ... Additional options passsed to \code{\link{ich_preprocess}}
+#'
+#' @return List of output registered and native space
+#' prediction/probability images
+#' @importFrom neurobase remake_img
+#' @importFrom extrantsr ants_bwlabel
+#' @import randomForest
+#' @seealso \code{\link{ich_candidate_voxels}}
+#' @export
+ich_predict = function(df,
+                       nim,
+                       model = c("rf", "logistic", "big_rf"),
+                       verbose = TRUE,
+                       native = TRUE,
+                       native_img = NULL,
+                       transformlist = NULL,
+                       interpolator = NULL,
+                       native_thresh = 0.5,
+                       shiny = FALSE,
+                       model_list = NULL,
+                       smoothed_cutoffs = NULL,
+                       outfile = NULL,
+                       ...) {
+
+  # if (!have_matlab()) {
+  #   stop("MATLAB Path not defined!")
+  # }
+
+  if (is.null(outfile)) {
+    outfile = tempfile(fileext = ".nii.gz")
+  }
+  cn = colnames(df)
+  if (!("multiplier" %in% cn)) {
+    df$multiplier = ich_candidate_voxels(df)
+  }
+  df$Y = NULL
+  cc = complete.cases(df)
+  if (!all(cc)) {
+    warning("NAs or missing in DF, removing")
+    for (icn in seq(ncol(df))) {
+      x = df[, icn]
+      if (!(class(x) %in% c("factor", "character"))) {
+        x[ !is.finite(x) ] = 0
+      }
+      df[, icn] = x
+    }
+  }
+  msg = "# Making Prediction"
+  if (verbose) {
+    message(msg)
+  }
+  if (shiny) {
+    shiny::incProgress(message = msg)
+  }
+  env = as.environment("package:ichseg")
+
+  # Getting modlist for model and cutoff
+  if (is.null(model_list)) {
+    modlist.name = paste0(model, "_modlist")
+    modlist = env[[modlist.name]]
+  } else {
+    modlist = model_list
+  }
+  mod = modlist$mod
+  cutoff = modlist$mod.dice.coef[1, "cutoff"]
+
+  rm(list = c("modlist"))
+
+  # Getting smoothed cutoff
+  if (is.null(smoothed_cutoffs)) {
+    smoothed_name = paste0("smoothed_", model, "_cutoffs")
+    scutoffs = env[[smoothed_name]]
+  } else {
+    scutoffs = smoothed_cutoffs
+  }
+
+  smoothed_cutoff = scutoffs$mod.dice.coef[1, "cutoff"]
+  rm(list = c("scutoffs", "smoothed_name"))
+
+  p = switch(model,
+             rf = predict(mod,
+                          newdata = df[ df$multiplier, ],
+                          type = "prob")[, "1"],
+             big_rf = predict(mod,
+                              newdata = df[ df$multiplier, ],
+                              type = "prob")[, "1"],
+             logistic = predict(mod,
+                                df[ df$multiplier, ],
+                                type = "response"))
+  msg = "# Making Prediction Image"
+  if (verbose) {
+    message(msg)
+  }
+  nim = check_nifti(nim)
+  mult_img = niftiarr(nim, df$multiplier)
+
+  # p = predict(mod, df[ df$multiplier, ], type = "response")
+  pimg = remake_img(p,
+                    nim,
+                    mult_img)
+
+  mask = niftiarr(nim, df$mask)
+  pimg = mask_img(pimg, mask)
+  msg = "# Smoothing Image"
+  if (verbose) {
+    message(msg)
+  }
+  sm.pimg  = mean_image(pimg,
+                        nvoxels = 1,
+                        verbose = verbose)
+  sm.pimg[abs(sm.pimg) <
+            .Machine$double.eps ^ 0.5 ] = 0
+  sm.pimg = niftiarr(nim, sm.pimg)
+  sm.pimg[is.na(sm.pimg)] = 0
+
+
+  sm.pred = sm.pimg > smoothed_cutoff
+  pred = pimg > cutoff
+
+  msg = "# Connected Components"
+  if (verbose) {
+    message(msg)
+  }
+  # cc = spm_bwlabel(pred, k = 100)
+  # scc = spm_bwlabel(sm.pred, k = 100)
+  cc = ants_bwlabel(img = pred, k = 100, binary = TRUE)
+  scc = ants_bwlabel(img = sm.pred, k = 100, binary = TRUE)
+
+  ##############################################################
+  # Back to Native Space!
+  ##############################################################
+  res = list(
+    prediction_image = cc,
+    smoothed_prediction_image = scc,
+    probability_image = pimg,
+    smoothed_probability_image = sm.pimg)
+
+  ##############################################################
+  # Inverted!
+  ##############################################################
+  native_res = NULL
+  if (native) {
+    msg = "# Projecting back to Native Space"
+    if (verbose) {
+      message(msg)
+    }
+    stopifnot(!is.null(interpolator))
+    stopifnot(!is.null(transformlist))
+    native_res = lapply(res, function(x){
+      ants_apply_transforms(fixed = native_img,
+                            moving = x,
+                            transformlist = transformlist,
+                            interpolator = interpolator,
+                            whichtoinvert = c(1)
+      )
+    })
+    native_res$smoothed_prediction_image = neurobase::datatyper(
+      native_res$smoothed_prediction_image > native_thresh
+    )
+    native_res$prediction_image = neurobase::datatyper(
+      native_res$prediction_image > native_thresh
+    )
+    writenii(native_res$smoothed_prediction_image, outfile)
+  } else {
+    writenii(res$smoothed_prediction_image, outfile)
+  }
+  res$cutoff = cutoff
+  res$smoothed_cutoff = smoothed_cutoff
+
+  L = list(registered_prediction = res)
+  L$native_prediction = native_res
+  L$outfile = outfile
+  return(L)
+}