R/maxent.predict.R

In correapvf/caretSDM: Species Distribution Models using caret

#' Methods for maxent object
#'
#' predict(), summary(), plot() and var_imp() methods for maxent object
#'
#' \itemize{
#' \item predict - return predicted classes or matrix of class probabilities
#' \item summary - open report generated by maxent.jar
#' \item plot - Plot a ROC curve based in the model
#' \item var_imp - Get variable importance for a maxent model
#' }
#' @param x,object A model returned by \code{\link{maxent}}.
#' @param newdata A data.frame of values to predict. If \code{NULL}, predictions of the training samples are returned.
#' @param thrtype Which thr type to use. Only used if \code{type = 'raw'}.
#' @param clamp logical. Apply clamping in prediction?
#' @param type One of 'raw' or 'prob'.
#' @param ... ignored
#' @name methods.maxent
#' @rdname methods.maxent
#' @export
predict.maxent <- function(object, newdata = NULL, thrtype = NULL, clamp = NULL, type = "raw", ...) {

    if (is.null(newdata)) {
        pred <- object$predicted

    } else {

        thereisfile <- TRUE
        while (thereisfile) {
            f <- paste(round(stats::runif(10) * 9), collapse = "")
            predfile <- paste0(object$filesPath, '/predIn', f)
            thereisfile <- file.exists(predfile)
        }

        outfile <- paste0(object$filesPath, '/predOut', f, '.asc')

        check_matrix <- is.matrix(newdata)
        newdata <- as.data.frame(newdata)[, object$xNames]

        # check factors
        if (!is.null(object$params$factors)) {
            if (check_matrix) newdata <- matrix.inverse(newdata, object$params$factors)
            fi <- colnames(newdata) %in% object$params$factors
            newdata[, fi] <- sapply(newdata[, fi], as.numeric)
            colnames(newdata)[fi] <- paste0("FACTOR_", object$params$factors)
            args0 <- c("-t", "FACTOR_", "writeclampgrid=false")
        } else {
            args0 <- "writeclampgrid=false"
        }

        # write pred file
        n <- nrow(newdata)
        newdata <- data.table::data.table(species = rep("sp", n),
                                          x = 1:n, y = 1:n, newdata)
        data.table::fwrite(newdata, file = predfile, quote = FALSE)

        # clamp
        if (is.null(clamp)) clamp <- object$params$clamp
        doclamp <- if (clamp) "doclamp=true" else "doclamp=false"



        # calculate predictions
        invisible(system2('java', args = c("-mx1024m", "-cp", object$maxentPath, "density.Project", object$lambda,
                                        predfile, outfile, paste0("outputformat=", tolower(object$params$outputType)),
                                        doclamp, args0)))

        pred <- data.table::fread(outfile, sep =  "\n", header = FALSE, skip = 6)[[1L]]
        if (nrow(newdata) > 1) pred <- as.numeric(strsplit(pred, " ")[[1]])
        unlink(c(predfile, outfile))
    }

    if (type == "prob") {
        result <- data.frame(pres = pred,
                             aus = 1 - pred)
        colnames(result) <- object$classes

    } else if (type == "raw") {
        if (is.null(thrtype)) {
            thrtype <- object$params$thr
        } else {
            thrtype <- maxentThr(thrtype)
        }

        thr <- object$results$value[which(object$results$variable == paste(thrtype, object$params$outputType, "threshold"))]
        result <- factor(pred >= thr, levels = c(TRUE, FALSE), labels = object$classes)

    } else {
        stop("Type must be 'prob' or 'raw'")
    }

    return(result)
}


#' @export
print.maxent <- function(x, ...) {
    cat("Object of class maxent\n")
    print(x$call)
    cat("Reg  :", x$params$reg.features, "\n")
    cat("Beta :", x$params$beta, "\n")
    cat("Path :", x$path, "\n")
}


#' @rdname methods.maxent
#' @export
summary.maxent <- function(object, ...) {
    utils::browseURL(paste0('file://', object$path, '/species.html'))
}



#' @param testdata A data.frame with test values to be use in ROC curve. If \code{NULL}, training values will be used.
#' @param testy A factor of response variable of \code{testdata}. If \code{NULL}, it will be guessed from
#' \code{testdata} data.frame.
#' @param ... Further arguments to \code{\link[pROC]{plot.roc}}.
#' @rdname methods.maxent
#' @export
plot.maxent <- function(x, testdata = NULL, testy = NULL, ...) {
    if (is.null(testdata)) {
        rocObject <- pROC::roc(x$y, x$predicted, levels = rev(x$classes),
                               direction = "auto", quiet = TRUE)
    } else {
        pred.name <- "testy"
        if (is.null(testy)) {
            if (!is.null(x$terms)) {
                pred.name <- all.vars(x$terms)[attr(x$terms, "response")]
            } else {
                stop("Either provide 'testy' or train using a formula")
            }
            testy <- testdata[[pred.name]]
        }
        if (!all(levels(testy) %in% x$classes))
            stop(pred.name, " must have levels ", paste(x$classes, collapse = " and "))

        rocObject <- pROC::roc(testy,
                               predict(x, testdata, type = "prob")[, x$classes[1]],
                               levels = rev(x$classes),
                               direction = "auto", quiet = TRUE)
    }
    return(invisible(pROC::plot.roc(rocObject, ...)))
}



#' @export
var_imp <- function(x, ...) UseMethod("var_imp")

#' @param itype One of 'contribution', 'permutation' or 'both'.
#' @rdname methods.maxent
#' @export
var_imp.maxent <- function(object, itype="both", ...) {
    index <- endsWith(object$results$variable, "contribution")
    if (itype == "contribution") {
        out <- data.frame(Overall = object$results$value[index])
    } else if (itype == "permutation") {
        out <- data.frame(perm.imp = object$results$value[endsWith(object$results$variable, "permutation importance")])
    } else {
        out <- data.frame(Overall = object$results$value[index],
                          perm.imp = object$results$value[endsWith(object$results$variable, "permutation importance")])
    }

    varnames <- sapply(strsplit(object$results$variable[index], " "), `[`, 1)
    rownames(out) <- sub("^FACTOR_", "", varnames)
    return(out)
}