V: Volume prediction model

In wslerry/tRee: Volume/Ha Prediction Model of Acacia Mg.

Usage

V(data, model)

Arguments

data
model

Examples

##---- Should be DIRECTLY executable !! ----
##-- ==>  Define data, use random,
##--	or do  help(data=index)  for the standard data sets.

## The function is currently defined as
function (data, model)
{
    library(randomForest)
    library(e1071)
    library(caret)
    library(ggplot2)
    ggplotRegression <- function(fit) {
        a <- signif(coef(fit)[1], digits = 5)
        b <- signif(coef(fit)[2], digits = 5)
        if (coef(fit)[2] >= 0) {
            textlab <- paste("y = ", a, " + ", b, "x", sep = "")
        }
        else {
            textlab <- paste("y = ", a, " - ", b, "x", sep = "")
        }
        options(repr.plot.width = 4, repr.plot.height = 4)
        ggplot(fit$model, aes_string(x = names(fit$model)[1],
            y = names(fit$model)[2])) + geom_point() + geom_smooth(method = "lm",
            col = "red", size = 0.5, se = TRUE) + labs(x = "Observations",
            y = "Predictions", title = paste("Adj. R2 = ", signif(summary(fit)$adj.r.squared,
                5), " | ", textlab)) + theme(plot.title = element_text(size = 8,
            face = "bold"))
    }
    traindata <- base::sample(nrow(data), size = 0.7 * nrow(data),
        replace = FALSE)
    TrainSet <- data[traindata, ]
    ValidSet <- data[-traindata, ]
    ctrl <- caret::trainControl(method = "cv", number = 5, savePredictions = TRUE)
    if (model == "mod1") {
        mod <- train(log(VUB) ~ I(log(DBH)) + I(log(H)) + I(log(BA)),
            data = TrainSet, method = "lm", trControl = ctrl,
            metric = "Rsquared")
        predictions <- predict(mod, ValidSet)
        predicted_V <- data.frame(log(ValidSet$VUB), predictions)
        fit.mod <- lm(log(ValidSet$VUB) ~ predictions, data = predicted_V)
        g <- ggplotRegression(fit.mod)
        print(mod$finalModel)
        g
    }
    else if (model == "mod2") {
        mod <- train(log(VUB) ~ I(log(DBH)) + I(log(H)) + I(log(DBH)^2) +
            I(log(H)^2), data = TrainSet, method = "lm", trControl = ctrl,
            metric = "Rsquared")
        predictions <- predict(mod, ValidSet)
        predicted_V <- data.frame(log(ValidSet$VUB), predictions)
        fit.mod <- lm(log(ValidSet$VUB) ~ predictions, data = predicted_V)
        g <- ggplotRegression(fit.mod)
        print(mod$finalModel)
        g
    }
    else if (model == "mod3") {
        mod <- train(log(VUB) ~ I(log(DBH^2 * H)), data = TrainSet,
            method = "lm", trControl = ctrl, metric = "Rsquared")
        predictions <- predict(mod, ValidSet)
        predicted_V <- data.frame(predictions, log(ValidSet$VUB))
        fit.mod <- lm(log(ValidSet$VUB) ~ predictions, data = predicted_V)
        g <- ggplotRegression(fit.mod)
        print(mod$finalModel)
        g
    }
    else if (model == "mod4") {
        set.seed(10)
        mod <- train(log(VUB) ~ I(log(DBH^2)) + I(log(H^2)),
            data = TrainSet, method = "lm", trControl = ctrl,
            metric = "Rsquared")
        predictions <- predict(mod, ValidSet)
        predicted_V <- data.frame(predictions, log(ValidSet$VUB))
        fit.mod <- lm(log(ValidSet$VUB) ~ predictions, data = predicted_V)
        g <- ggplotRegression(fit.mod)
        print(mod$finalModel)
        g
    }
    else if (model == "mod5") {
        set.seed(10)
        mod <- train(log(VUB) ~ I(log(H)), data = TrainSet, method = "lm",
            trControl = ctrl, metric = "Rsquared")
        predictions <- predict(mod, ValidSet)
        predicted_V <- data.frame(predictions, log(ValidSet$VUB))
        fit.mod <- lm(log(ValidSet$VUB) ~ predictions, data = predicted_V)
        g <- ggplotRegression(fit.mod)
        print(mod$finalModel)
        g
    }
    else if (model == "modRF") {
        set.seed(10)
        mod <- train(log(VUB) ~ I(log(H)), data = TrainSet, method = "rf",
            ntree = 1000, trControl = ctrl, metric = "Rsquared",
            importance = T, proximity = T)
        predictions <- predict(mod, ValidSet)
        predicted <- data.frame(predictions, log(ValidSet$VUB))
        fit.mod <- lm(log(ValidSet$VUB) ~ predictions, data = predicted)
        g <- ggplotRegression(fit.mod)
        fit.mod$coefficients
        print(mod$finalModel)
        g
    }
    else if (model == "modRF2") {
        rf <- randomForest(VUB ~ H, TrainSet, ntree = 1000, importance = T,
            proximity = TRUE, keep.forest = TRUE)
        print(rf)
        plot(rf, main = "Random Forest")
        predictions <- predict(rf, ValidSet)
        predicted <- data.frame(log(predictions), log(ValidSet$VUB))
        fit.mod <- lm(log(ValidSet$VUB) ~ log(predictions), data = predicted)
        plot(predicted) + abline(fit.mod)
    }
    else {
        stop("No Model selected. Please select model eg., 'mod1','mod2','mod3','mod4','mod5' or 'modRF'")
    }
  }

wslerry/tRee documentation built on Oct. 5, 2019, 10:30 p.m.