---
title: "Predicting Pine Beetle Brood Tree Distance From the Pine Beetle Data Set"
output:
  flexdashboard::flex_dashboard:
    orientation: rows
    vertical_layout: scroll
    theme: 
      version: 4
      bootswatch: minty
    source_code: embed
    #https://rstudio.github.io/flexdashboard/articles/theme.html
---

```{r setup, include=FALSE}
library(flexdashboard)
library(tidymodels)
library(tidyverse)
library(vip)
library(plotly)
library(gt)
library(readxl)
library(bslib)
library(ranger)
theme_set(theme_minimal())

pine_beetle <- read_xlsx(
  "~/Desktop/Spring Semester 2026/DataScience2/hgen-612/data/Data_1993.xlsx")
#pine_beetle <- read_xlsx("data/Data_1993.xlsx") # For reproducibility purposes
```

```{r data-prep}
pine_data <- pine_beetle %>%
  drop_na(DeadDist)

set.seed(123)

split <- initial_split(pine_data, prop = 0.8)

train <- training(split)
test  <- testing(split)
```

```{r recipe}
# make a recipe
pine_recipe <- recipe(DeadDist ~ ., data = train) %>%
  step_zv(all_predictors()) %>% # remove vars with no variation
  step_impute_mean(all_numeric_predictors()) %>%
  step_normalize(all_numeric_predictors())
```

```{r models}
lm_model <- linear_reg() %>%
  set_engine("lm")

rf_model <- rand_forest( #using a random forest model due to complexity of data, IDK much about ecology and I expect the var relationships to not all be linear.
  trees = 500,
  min_n = 5
) %>%
  set_engine("ranger") %>%
  set_mode("regression")
```

```{r workflows}
# compare lm and rf model
lm_workflow <- workflow() %>%
  add_recipe(pine_recipe) %>%
  add_model(lm_model)

rf_workflow <- workflow() %>%
  add_recipe(pine_recipe) %>%
  add_model(rf_model)

lm_fit <- fit(lm_workflow, train)
rf_fit <- fit(rf_workflow, train)
```

```{r predictions}
lm_preds <- predict(lm_fit, test) %>%
  bind_cols(test)

rf_preds <- predict(rf_fit, test) %>%
  bind_cols(test)
```

Page 1
=======================================================================

Column
-----------------------------------------------------------------------

### Dataset Overview
*Sample of predictors used to generate the brood tree distance model.*

```{r}
pine_data %>%
  sample_n(10) %>%
  gt()
```

Column
-----------------------------------------------------------------------

### Distribution of Distance to Brood Tree
*Understanding the response variable.*

```{r}
dist_plot <- pine_data %>%
  ggplot(aes(x = DeadDist)) + 
  geom_histogram(bins = 30, fill = "#2C7FB8", alpha = 0.8) +
  labs(
    title = "Distribution of Distance to Nearest Brood Tree",
    x = "Distance (DeadDist)",
    y = "Count"
  )

ggplotly(dist_plot) %>%
  layout(margin = list(l = 70, r = 30, t = 60, b = 100)) # created bigger margins b/c redering was cutting off my axis
```

Column
-----------------------------------------------------------------------

### Tree Diameter vs Distance
*Larger trees potentially influence brood tree spacing.*

```{r}
diam_plot <- pine_data %>%
  ggplot(aes(TreeDiam, DeadDist)) +
  geom_point(alpha = 0.4) +
  geom_smooth(method = "lm") +
  labs(
    title = "Tree Diameter vs Distance to Brood Tree",
    x = "Diameter of Tree",
    y = "Distance to Brood Tree"
  )

ggplotly(diam_plot) %>%
  ggplotly(dist_plot) %>%
  layout(margin = list(l = 70, r = 30, t = 60, b = 100))
```

Page 2
=======================================================================

Column
-----------------------------------------------------------------------

### Linear Regression Model Results
*Predicted vs observed distances.*

```{r}
lm_plot <- lm_preds %>%
  ggplot(aes(DeadDist, .pred)) +
  geom_point(alpha = 0.4) +
  geom_abline(linetype = 2, color = "red") +
  labs(
    title = "Linear Regression Predictions",
    x = "Observed Distance",
    y = "Predicted Distance"
  )

ggplotly(lm_plot) %>%
  ggplotly(dist_plot) %>%
  layout(margin = list(l = 70, r = 30, t = 60, b = 100))
```

Column
-----------------------------------------------------------------------

### Random Forest Model Results
*Comparing nonlinear prediction and its performance.*

```{r}
rf_plot <- rf_preds %>%
  ggplot(aes(DeadDist, .pred)) +
  geom_point(alpha = 0.4) +
  geom_abline(linetype = 2, color = "red") +
  labs(
    title = "Random Forest Predictions",
    x = "Observed Distance",
    y = "Predicted Distance"
  )

ggplotly(rf_plot) %>%
  ggplotly(dist_plot) %>%
  layout(margin = list(l = 70, r = 30, t = 60, b = 100))
```

Column
-----------------------------------------------------------------------

### Comparing Performance of the Models 
*Random forest model improves predictive accuracy.*

```{r}
lm_metrics <- lm_preds %>%
  metrics(truth = DeadDist, estimate = .pred) %>%
  mutate(model = "Linear Regression")

rf_metrics <- rf_preds %>%
  metrics(truth = DeadDist, estimate = .pred) %>%
  mutate(model = "Random Forest")

bind_rows(lm_metrics, rf_metrics) %>%
  select(model, .metric, .estimate) %>%
  pivot_wider(names_from = .metric, values_from = .estimate) %>%
  gt()
```

Page 3
=======================================================================

Column
-----------------------------------------------------------------------

### Residual Diagnostics
*Evaluating model errors.*

```{r}
res_plot <- lm_preds %>%
  mutate(residual = DeadDist - .pred) %>%
  ggplot(aes(.pred, residual)) +
  geom_point(alpha = 0.4) +
  geom_hline(yintercept = 0, linetype = 2) +
  labs(
    title = "Residual Diagnostics",
    x = "Predicted Distance",
    y = "Residual"
  )

ggplotly(res_plot) %>%
  layout(margin = list(l = 70, r = 30, t = 60, b = 100))
```

Column
-----------------------------------------------------------------------

### QQ Plot of Residuals
*Assessing normality.*

``` {r}
qq_plot <- lm_preds %>%
  mutate(residual = DeadDist - .pred) %>%
  ggplot(aes(sample = residual)) +
  stat_qq(alpha = 0.5) +
  stat_qq_line(color = "red", linetype = 2) +
  labs(
    title = "QQ Plot of Residuals",
    x = "Theoretical Quantiles",
    y = "Sample Quantiles"
  )

qq_plot 
```