I've built a shiny app to model 3D data using Thin Plate Splines, and Tensor Product Smooths.
Note: Both models are prefit to speed up plot selection (input$Mod
)
Unfortunately, when I try call input$Mod
using functions such as get()
it breaks this line;
Z <- matrix(predict(Model2b(), newdat), steps, steps)
Example
Z <- matrix(predict(get(input$Mod), newdat), steps, steps)
How can I call model fits? I don't want to remodel the same data on repeat each time the plot choice is changed.
The app reads a locally stored CSV
Get data
mtcars <- get(data("mtcars"))
write.csv(mtcars, file = 'mtcars".csv', row.names = FALSE)
Inside the app please only select continuous numerical values such as hp, mpg, qsec, wt, disp, etc. Otherwise the splines won't fit.
Shiny App
# Clear all
rm(list = ls(all.names = T))
gc()
library(tidyverse)
library(dlookr)
library(scales)
library(shinyWidgets)
library(rgl)
library(mgcv)
# UI
ui <- navbarPage(title = "Advanced Analytics",
tabPanel("GAM Models",
sidebarLayout(
sidebarPanel(width = 3,
fileInput(inputId = "file1",
label = "Upload CSV",
accept = c(".csv")),
uiOutput("RespSelector"),
uiOutput("PredSelector"),
sliderTextInput(
inputId = "OutlierStrength",
label = "Outlier removal strength",
selected = 1,
choices = c('None',1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20)
),
numericInput(inputId = "ChangeK_TPS",
label = "Flexibility: Thin Plate Spline",
value = 0,
min = -50,
max = 50,
step = 1),
numericInput(inputId = "ChangeK_Te",
label = "Flexibility: Tensor Product Smooth",
value = 0,
min = -50,
max = 50,
step = 1),
verbatimTextOutput("aic"),
selectInput(inputId = "Mod",
label = "Model Type:",
choices = c("Thin Plate Spline" = 'Model1b()',
"Tensor Product Smooth" = 'Model2b()'))
),
mainPanel(
rglwidgetOutput("myplot", width = "1000px", height = "1280px"),
# DT::dataTableOutput(outputId = "Table1", width = "125%"),
br(),
br(),
br(),
br(),
br(),
verbatimTextOutput("Summary1"),
verbatimTextOutput("Summary2")
)
)
)
)
# Server
server <- function(input, output, session) {
# Upload CSV
csv_data <- reactive({req(input$file1)
# Read CSV and lightly trim tails
read_csv(input$file1$datapath) %>%
rowid_to_column("ID")
})
# Extract numeric colnames
VARS_numeric <- reactive({req(csv_data())
csv_data() %>%
select(where(is.numeric), -ID) %>%
colnames() %>%
sort()
})
# Render respone for UI selector
output$RespSelector <- renderUI({req(csv_data())
selectizeInput(inputId = "response",
label = "Select 1 response variable",
selected = NULL,
choices = VARS_numeric(),
multiple = FALSE)
})
# Render predictor UI selector
output$PredSelector <- renderUI({req(csv_data())
selectizeInput(inputId = "predictors",
label = "Select 2 predictors variables",
choices = VARS_numeric()[!(VARS_numeric() %in% input$response)],
multiple = TRUE,
options = list(maxItems = 2))
})
# Pivot
Pivot <- reactive({req(input$predictors, input$response, VARS_numeric(), csv_data())
csv_data() %>%
pivot_longer(cols = c(input$predictors, input$response)) %>%
arrange(ID)
})
# Normality test
NormalityTest <- reactive({req(Pivot())
Pivot() %>%
group_by(name) %>%
normality(value) %>%
ungroup() %>%
mutate(NormalDist = if_else(p_value <= 0.10, 'No', 'Yes')) %>%
select(name, NormalDist) %>%
left_join(., Pivot(), by = c('name' = 'name'))
})
# Find distribution parameters and create outlier limits
DistParameters <- reactive({req(NormalityTest())
NormalityTest() %>%
drop_na() %>%
group_by(name) %>%
summarise(across(value, list(
Median = ~median(value, na.rm = T),
Mean = ~mean(value, na.rm = T, trim = 0.05),
SD = ~sd(value, na.rm = T),
MAD = ~sum(abs(value - median(value, na.rm = T)))/length(value)))) %>%
ungroup() %>%
rename_with(~str_remove(., 'value_')) %>%
mutate(LowerLimit_MAD = Median - (MAD * ifelse(input$OutlierStrength == 'None', 99999, (as.numeric(input$OutlierStrength) - 21) *-1 /8 2) )) %>%
mutate(UpperLimit_MAD = Median (MAD * ifelse(input$OutlierStrength == 'None', 99999, (as.numeric(input$OutlierStrength) - 21) *-1 /8 2) )) %>%
mutate(LowerLimit_SD = Mean - (SD * ifelse(input$OutlierStrength == 'None', 99999, (as.numeric(input$OutlierStrength) - 21) *-1 /8 2) )) %>%
mutate(UpperLimit_SD = Mean (SD * ifelse(input$OutlierStrength == 'None', 99999, (as.numeric(input$OutlierStrength) - 21) *-1 /8 2) ))
})
# Remove outliers
NoOutliers <- reactive({req(DistParameters(), NormalityTest())
inner_join(NormalityTest(), DistParameters()) %>%
mutate(OutlierMAD = ifelse(NormalDist == 'No' & (value < LowerLimit_MAD | value > UpperLimit_MAD), 'Yes', 'No')) %>%
mutate(OutlierSD = ifelse(NormalDist == 'Yes' & (value < LowerLimit_SD | value > UpperLimit_SD), 'Yes', 'No')) %>%
mutate(Outlier = ifelse(NormalDist == 'No', OutlierMAD, OutlierSD)) %>% # Combine with row below
filter(Outlier == 'No') %>%
select(ID, name, value) %>%
pivot_wider(names_from = name, values_from = c(value)) %>%
select(-ID) %>%
drop_na()
})
# s(x1,x2) modela
ModelEquation1a <- reactive({req(NoOutliers(), input$response, input$predictors[1], input$predictors[2])
Equation1a <- as.formula(paste0(input$response," ~ ", 's(', input$predictors[1],',', input$predictors[2], ', bs = "tp")'))
})
# te(x1,x2) modela
ModelEquation2a <- reactive({req(NoOutliers(), input$response, input$predictors[1], input$predictors[2])
Equation2a <- as.formula(paste0(input$response,' ~ ', 'te(',input$predictors[1],',',input$predictors[2],')'))
})
# Model 1a
Model1a <- reactive({req(ModelEquation1a(), NoOutliers(), input$response, input$predictors[1], input$predictors[2])
gam(ModelEquation1a(), method="REML", data = NoOutliers())
})
# Model 2a
Model2a <- reactive({req(ModelEquation2a(), NoOutliers(), input$response, input$predictors[1], input$predictors[2])
gam(ModelEquation2a(), method="REML", data = NoOutliers())
})
# s(x1,x2) modelb
ModelEquation1b <- reactive({req(Model1a(), NoOutliers(), input$response, input$predictors[1], input$predictors[2], input$ChangeK_TPS)
Equation1b <- as.formula(paste0(input$response," ~ ", 's(', input$predictors[1],',', input$predictors[2], ', bs = "tp", k = ', Model1a()$rank input$ChangeK_TPS, ')'))
})
# te(x1,x2) modelb
ModelEquation2b <- reactive({req(Model2a(), NoOutliers(), input$response, input$predictors[1], input$predictors[2], input$ChangeK_Te)
Equation2b <- as.formula(paste0(input$response,' ~ ', 'te(',input$predictors[1],',',input$predictors[2],', k = ', sqrt(Model2a()$rank) input$ChangeK_Te, ')'))
})
# Model 1b
Model1b <- reactive({req(ModelEquation1b())
gam(ModelEquation1b(), method="REML", data = NoOutliers())
})
# Model 2b
Model2b <- reactive({req(ModelEquation2b())
gam(ModelEquation2b(), method="REML", data = NoOutliers())
})
# Summary1
output$Summary1 <- renderPrint({req(Model1b())
summary(Model1b())
})
# Summary2
output$Summary2 <- renderPrint({req(Model2b())
summary(Model2b())
})
# AIC
output$aic <- renderPrint({req(Model1b(), Model2b())
AIC(Model1b(), Model2b()) %>%
rownames_to_column(var = "Model Type") %>%
mutate(`Model Type` = replace(`Model Type`, str_detect(`Model Type`, "Model1b()"), "Thin Plate Spline")) %>%
mutate(`Model Type` = replace(`Model Type`, str_detect(`Model Type`, "Model2b()"), "Tensor Product Smooth")) %>%
select(`Model Type`, AIC) %>%
arrange(AIC)
})
# Data table Tab-1
output$Table1 <- DT::renderDataTable({
DT::datatable(data = NoOutliers(),
options = list(pageLength = 20),
rownames = FALSE)
})
output$myplot <- renderRglwidget({req(Model1b(), Model2b(), input$Mod, NoOutliers(), input$response, input$predictors[1], input$predictors[2], input$ChangeK_Te, input$ChangeK_TPS, input$OutlierStrength)
# Expand to grid
steps <- 30
FlyDF <- NoOutliers()
X1 <- with(FlyDF, seq(min(FlyDF[input$predictors[1]], na.rm = T), max(FlyDF[input$predictors[1]], na.rm = T), length = steps))
X2 <- with(FlyDF, seq(min(FlyDF[input$predictors[2]], na.rm = T), max(FlyDF[input$predictors[2]], na.rm = T), length = steps))
newdat <- expand.grid(X1 = X1, X2 = X2)
colnames(newdat) <- c(as.character(input$predictors[1]), as.character(input$predictors[2]))
Z <- matrix(predict(Model2b(), newdat), steps, steps)
# Plot
try(close3d())
persp3d(X1, X2, Z, col="yellow", polygon_offset = 1, axes = TRUE, box = TRUE,
xlab = as.character(input$predictors[1]), ylab = as.character(input$predictors[2]), zlab = as.character(input$response))
persp3d(X1, X2, Z, front = "lines", back = "lines", lit = FALSE, add = TRUE)
points3d(FlyDF[,c(input$predictors[1], input$predictors[2], input$response)], col = "blue")
rglwidget()
})
}
# Create Shiny app
shinyApp(ui = ui, server = server)
CodePudding user response:
If that is the only issue, you change your selectInput()
as
selectInput(inputId = "Mod", label = "Model Type:", choices = c("Thin Plate Spline" = 'Model1',
"Tensor Product Smooth" = 'Model2'))
Then create a eventReactive model as
myModel <- eventReactive(input$Mod, {
switch(input$Mod,
"Model1" = Model1b(),
"Model2" = Model2b())
})
and lastly use this in predict as
Z <- matrix(predict(myModel(), newdat), steps, steps)