inst/doc/getting-started-with-cheem.R

In cheem: Interactively Explore Local Explanations with the Radial Tour

## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
  echo       = TRUE,   # code
  include    = TRUE,   # plots
  results    = "show", # text: "hide", "show"
  ## stop("REPLACE ME"):
  eval       = FALSE,   # chunk code
  message    = FALSE,
  warning    = FALSE,
  error      = FALSE,
  collapse   = TRUE,
  comment    = "#>",
  fig.height = 4,
  fig.width  = 6,
  fig.align  = "center",
  cache      = FALSE
)

## ----limenonlinear, echo=FALSE, fig.cap="Illustration of non-linear classification boundary. The use of local explanations approximates the feature importance in the vicinity of one instance. This allow us to understand a change in which features would result in a red plus being classified as a blue circle. From _Ribiro, M. et. all. (2017). Why should I trust you?_"----
# knitr::include_graphics("../inst/shiny_apps/cheem/www/lime_nonlinear.png")

## -----------------------------------------------------------------------------
# ## Download if not installed
# if(!require(cheem))    install.packages("cheem", dependencies = TRUE)
# if(!require(treeshap)) install.packages("treeshap", dependencies = TRUE)
# if(!require(shapviz))  install.packages("shapviz", dependencies = TRUE)
# ## Load onto session
# library(cheem)
# library(xgboost)
# library(shapviz)
# 
# ## Setup
# X    <- amesHousing2018_NorthAmes[, 1:9]
# Y    <- amesHousing2018_NorthAmes$SalePrice
# clas <- amesHousing2018_NorthAmes$SubclassMS
# 
# ## Model and predict
# ames_train    <- data.matrix(X) %>% xgb.DMatrix(label = Y)
# ames_xgb_fit  <- xgboost(data = ames_train, max.depth = 3, nrounds = 25)
# ames_xgb_pred <- predict(ames_xgb_fit, newdata = ames_train)
# ames_xgb_pred %>% head()
# 
# ## SHAP values
# shp <- shapviz(ames_xgb_fit, X_pred = ames_train, X = X)
# ## Keep just the [n, p] local explanations
# ames_xgb_shap <- shp$S
# ames_xgb_shap %>% head()

## -----------------------------------------------------------------------------
# ## Preprocessing for cheem analysis
# ames_chm <- cheem_ls(X, Y,
#                      class      = clas,
#                      attr_df    = ames_xgb_shap,
#                      pred       = ames_xgb_pred,
#                      label      = "Ames, xgb, shap")
# names(ames_chm)

## ----out.width="100%"---------------------------------------------------------
# prim <- 1
# comp <- 17
# global_view(ames_chm, primary_obs = prim, comparison_obs = comp,
#             height_px = 240, width_px = 720,
#             as_ggplot = TRUE, color = "log_maha.data")

## ----out.width="100%", eval = FALSE-------------------------------------------
# ## Normalized attribution basis of the PI
# bas <- sug_basis(ames_xgb_shap, rownum = prim)
# ## Default feature to manipulate:
# #### the feature with largest separation between PI and CI attribution
# mv  <- sug_manip_var(
#   ames_xgb_shap, primary_obs = prim, comparison_obs = comp)
# ## Make the radial tour
# ggt <- radial_cheem_tour(
#   ames_chm, basis = bas, manip_var = mv,
#   primary_obs = prim, comparison_obs = comp, angle = .15)
# 
# ## Animate it
# animate_gganimate(ggt, fps = 6)
#   #height = 2, width = 4.5, units = "in", res = 150
# ## Or as a plotly html widget
# #animate_plotly(ggt, fps = 6)

## ----echo=FALSE, out.width="100%"---------------------------------------------
# ## To mitigate file size (CRAN note) and run time create a gif and include that instead of executing code to make inline.
# if(FALSE){
#   prim <- 1
#   comp <- 17
#   bas <- sug_basis(shap_df, rownum = prim)
#   mv  <- sug_manip_var(
#     shap_df, primary_obs = prim, comparison_obs = comp)
#   ggt <- radial_cheem_tour(
#     this_ls, basis = bas, manip_var = mv,
#     primary_obs = prim, comparison_obs = comp, angle = .15)
# 
#   #### .gif is about .2 Mb saved, while HTML widget was about 7 Mb.
#   anim <- animate_gganimate(
#     ggt, fps = 6, height = 2, width = 4.5, units = "in", res = 150)
#   gganimate::anim_save("tour_penguins.gif", animation = anim)#, path = "./vignettes")
#   beepr::beep()
# }
# #knitr::include_graphics("tour_penguins.gif")

## -----------------------------------------------------------------------------
# library(ggplot2)
# prim <- 1
# 
# ggplot(penguins_na.rm, aes(x = bill_length_mm,
#                            y = flipper_length_mm,
#                            colour = species,
#                            shape = species)) +
#   geom_point() +
#   ## Highlight PI, *
#   geom_point(data = penguins_na.rm[prim, ],
#              shape = 8, size = 5, alpha = 0.8) +
#   ## Theme, scaling, color, and labels
#   theme_bw() +
#   theme(aspect.ratio = 1) +
#   scale_color_brewer(palette = "Dark2") +
#   labs(y = "Flipper length [mm]", x = "Bill length [mm]",
#        color = "Observed species", shape = "Observed species")

## ----eval=FALSE, echo=TRUE----------------------------------------------------
# if(!require(shapviz)) install.packages("shapviz")
# if(!require(xgboost)) install.packages("xgboost")
# library(shapviz)
# library(xgboost)
# set.seed(3653)
# 
# ## Setup
# X    <- spinifex::penguins_na.rm[, 1:4]
# Y    <- spinifex::penguins_na.rm$species
# clas <- spinifex::penguins_na.rm$species
# 
# ## Model and predict
# peng_train    <- data.matrix(X) %>%
#   xgb.DMatrix(label = Y)
# peng_xgb_fit  <- xgboost(data = peng_train, max.depth = 3, nrounds = 25)
# peng_xgb_pred <- predict(peng_xgb_fit, newdata = peng_train)
# 
# ## SHAP
# peng_xgb_shap <- shapviz(peng_xgb_fit, X_pred = peng_train, X = X)
# ## Keep just the [n, p] local explanations
# peng_xgb_shap <- peng_xgb_shap$S

## ----eval=FALSE, echo=TRUE----------------------------------------------------
# if(!require(treeshap)) install.packages("treeshap")
# if(!require(randomForest)) install.packages("randomForest")
# library(treeshap)
# library(randomForest)
# 
# ## Setup
# X    <- spinifex::wine[, -1:2]
# Y    <- spinifex::wine$Alcohol
# clas <- spinifex::wine$Type
# 
# ## Fit randomForest::randomForest
# wine_rf_fit  <- randomForest::randomForest(
#   X, Y, ntree = 125,
#   mtry = ifelse(is_discrete(Y), sqrt(ncol(X)), ncol(X) / 3),
#   nodesize = max(ifelse(is_discrete(Y), 1, 5), nrow(X) / 500))
# wine_rf_pred <- predict(wine_rf_fit)
# 
# ## treeshap::treeshap()
# wine_rf_tshap <- wine_rf_fit %>%
#   treeshap::randomForest.unify(X) %>%
#   treeshap::treeshap(X, interactions = FALSE, verbose = FALSE)
# ## Keep just the [n, p] local explanations
# wine_rf_tshap <- wine_rf_tshap$shaps

## ----eval=FALSE, echo=TRUE----------------------------------------------------
# if(!require(DALEX)) install.packages("DALEX")
# library(DALEX)
# 
# ## Setup
# X    <- dragons[, c(1:4, 6)]
# Y    <- dragons$life_length
# clas <- dragons$colour
# 
# ## Model and predict
# drag_lm_fit  <- lm(data = data.frame(Y, X), Y ~ .)
# drag_lm_pred <- predict(drag_lm_fit)
# 
# ## SHAP via DALEX, versatile but slow
# drag_lm_exp <- explain(drag_lm_fit, data = X, y = Y,
#                        label = "Dragons, LM, SHAP")
# ## DALEX::predict_parts_shap is flexible, but slow and one row at a time
# drag_lm_shap <- matrix(NA, nrow(X), ncol(X))
# sapply(1:nrow(X), function(i){
#   pps <- predict_parts_shap(drag_lm_exp, new_observation = X[i, ])
#   ## Keep just the [n, p] local explanations
#   drag_lm_shap[i, ] <<- tapply(
#     pps$contribution, pps$variable, mean, na.rm = TRUE) %>% as.vector()
# })
# drag_lm_shap <- as.data.frame(drag_lm_shap)