+Number of missing data per variable, shows that bp.2d and bp.2s have more than 50% missing entries
+
+
+
+
+
+
+
+Heatmap visualizing Pearson correlation coefficient between numerical variables
+
+
+
+
+
+
+
2.2 Data splitting
+
+
# use tidymodels framework to fit Lasso regression model for predicting BMI
+# using repeated cross-validation to tune lambda value in L1 penalty term
+
+# select random seed value
+myseed <-123
+
+# split data into non-test (other) and test (80% s)
+set.seed(myseed)
+data_split <-initial_split(data_diabetes_narm, strata = BMI, prop =0.8) # holds splitting info
+data_other <- data_split %>%training() # creates non-test set (function is called training but it refers to non-test part)
+data_test <- data_split %>%testing() # creates test set
+
+# prepare repeated cross-validation splits with 5 folds repeated 3 times
+set.seed(myseed)
+data_folds <-vfold_cv(data_other,
+v =5,
+repeats =3,
+strata = BMI)
+
+# check the split
+dim(data_diabetes)
+## [1] 403 20
+dim(data_other)
+## [1] 291 18
+dim(data_test)
+## [1] 75 18
+
+# check BMI distributions in data splits
+par(mfrow=c(3,1))
+hist(data_diabetes$BMI, xlab ="", main ="BMI: all", 50)
+hist(data_other$BMI, xlab ="", main ="BMI: non-test", 50)
+hist(data_test$BMI, xlab ="", main ="BMI: test", 50)
+
+
+
+
+Distribution of BMI values given all data and spits into non-test and test
+
+
+
+
+
+
+
2.3 Feature engineering
+
+
# create data recipe (feature engineering)
+
+inch2m <-2.54/100
+pound2kg <-0.45
+
+data_recipe <-recipe(BMI ~ ., data = data_other) %>%
+update_role(id, new_role ="sampleID") %>%
+step_mutate(height = height * inch2m, height =round(height, 2)) %>%# convert height to meters
+step_mutate(weight = weight * pound2kg, weight =round(weight, 2)) %>%# convert weight to kg
+step_rename(glu = stab.glu) %>%# rename stab.glu to glu
+step_log(glu) %>%#ln transform glucose
+step_zv(all_numeric()) %>%# removes variables that are highly sparse and unbalanced (if found)
+step_corr(all_numeric(), -all_outcomes(), -has_role("sampleID"), threshold =0.8) %>%# removes variables with large absolute correlations with other variables (if found)
+step_dummy(location, gender, frame) %>%# convert categorical variables to dummy variables
+step_normalize(all_numeric(), -all_outcomes(), -has_role("sampleID"), skip =FALSE)
+
+# you can implement more steps: see https://recipes.tidymodels.org/reference/index.html
+
+# print recipe
+data_recipe
+
+# check if recipe is doing what it is supposed to do
+# i.e. bake the data
+data_other_prep <- data_recipe %>%
+prep() %>%
+bake(new_data =NULL)
+
+## bake test data
+data_test_prep <- data_recipe %>%
+prep() %>%
+bake(new_data = data_test)
+
+# preview baked data
+print(head(data_other_prep))
+## # A tibble: 6 × 17
+## id chol glu hdl ratio glyhb age height bp.1s bp.1d hip
+## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
+## 1 1045 -0.319 -0.539 -0.830 0.486 -0.172 -0.645 -0.474 -1.16 -0.552 -1.66
+## 2 1271 0.449 -1.09 -0.324 0.320 -0.458 -1.39 -1.29 -1.59 -1.00 -0.914
+## 3 1277 -0.657 -0.572 2.32 -1.45 -0.642 -0.337 1.56 0.286 2.15 -1.29
+## 4 1303 -0.590 -0.410 -0.436 -0.178 -0.518 0.341 0.545 -0.309 0.500 -1.47
+## 5 1309 -0.206 -0.347 0.687 -0.730 -0.860 -1.32 0.0357 -0.735 -0.402 -1.66
+## 6 1315 -0.793 -0.572 0.350 -0.841 0.225 0.649 1.26 -0.565 -1.45 -1.29
+## # ℹ 6 more variables: time.ppn <dbl>, BMI <dbl>, location_Louisa <dbl>,
+## # gender_male <dbl>, frame_medium <dbl>, frame_small <dbl>
+
+
+
+
2.4 Lasso regression
+
+
# define model
+model <-linear_reg(penalty =tune(), mixture =1) %>%
+set_engine("glmnet") %>%
+set_mode("regression")
+
+# create workflow with data recipe and model
+wf <-workflow() %>%
+add_model(model) %>%
+add_recipe(data_recipe)
+
+# define parameters range for tuning
+grid_lambda <-grid_regular(penalty(), levels =25)
+
+# tune lambda
+model_tune <- wf %>%
+tune_grid(resamples = data_folds,
+grid = grid_lambda)
+
+# show metrics average across folds
+model_tune %>%
+collect_metrics(summarize =TRUE)
+## # A tibble: 50 × 7
+## penalty .metric .estimator mean n std_err .config
+## <dbl> <chr> <chr> <dbl> <int> <dbl> <chr>
+## 1 1 e-10 rmse standard 2.48 15 0.142 Preprocessor1_Model01
+## 2 1 e-10 rsq standard 0.851 15 0.0143 Preprocessor1_Model01
+## 3 2.61e-10 rmse standard 2.48 15 0.142 Preprocessor1_Model02
+## 4 2.61e-10 rsq standard 0.851 15 0.0143 Preprocessor1_Model02
+## 5 6.81e-10 rmse standard 2.48 15 0.142 Preprocessor1_Model03
+## 6 6.81e-10 rsq standard 0.851 15 0.0143 Preprocessor1_Model03
+## 7 1.78e- 9 rmse standard 2.48 15 0.142 Preprocessor1_Model04
+## 8 1.78e- 9 rsq standard 0.851 15 0.0143 Preprocessor1_Model04
+## 9 4.64e- 9 rmse standard 2.48 15 0.142 Preprocessor1_Model05
+## 10 4.64e- 9 rsq standard 0.851 15 0.0143 Preprocessor1_Model05
+## # ℹ 40 more rows
+
+# plot k-folds results across lambda range
+model_tune %>%
+collect_metrics() %>%
+ dplyr::filter(.metric =="rmse") %>%
+ggplot(aes(penalty, mean, color = .metric)) +
+geom_errorbar(aes( ymin = mean - std_err, ymax = mean + std_err), alpha =0.5) +
+scale_x_log10() +
+geom_line(linewidth =1.5) +
+theme_bw() +
+theme(legend.position ="none") +
+scale_color_brewer(palette ="Set1")
+
+# best lambda value (min. RMSE)
+model_best <- model_tune %>%
+select_best(metric ="rmse")
+
+print(model_best)
+## # A tibble: 1 × 2
+## penalty .config
+## <dbl> <chr>
+## 1 0.0562 Preprocessor1_Model22
+
+# finalize workflow with tuned model
+wf_final <- wf %>%
+finalize_workflow(model_best)
+
+# last fit
+fit_final <- wf_final %>%
+last_fit(split = data_split)
+
+# final predictions
+y_test_pred <- fit_final %>%collect_predictions() # predicted BMI
+
+# final predictions: performance on test (unseen data)
+fit_final %>%collect_metrics()
+## # A tibble: 2 × 4
+## .metric .estimator .estimate .config
+## <chr> <chr> <dbl> <chr>
+## 1 rmse standard 2.79 Preprocessor1_Model1
+## 2 rsq standard 0.857 Preprocessor1_Model1
+
+# plot predictions vs. actual for test data
+plot(data_test$BMI, y_test_pred$.pred, xlab="BMI (actual)", ylab ="BMI (predicted)", las =1, pch =19)
+
+# correlation between predicted and actual BMI values for test data
+cor(data_test$BMI, y_test_pred$.pred)
+## [1] 0.9256857
+
+# re-fit model on all non-test data
+model_final <- wf_final %>%
+fit(data_other)
+
+# show final model
+tidy(model_final)
+## # A tibble: 16 × 3
+## term estimate penalty
+## <chr> <dbl> <dbl>
+## 1 (Intercept) 28.7 0.0562
+## 2 chol 0 0.0562
+## 3 glu 0.0229 0.0562
+## 4 hdl 0 0.0562
+## 5 ratio 0.335 0.0562
+## 6 glyhb -0.0512 0.0562
+## 7 age -0.257 0.0562
+## 8 height -1.86 0.0562
+## 9 bp.1s -0.294 0.0562
+## 10 bp.1d 0.203 0.0562
+## 11 hip 5.60 0.0562
+## 12 time.ppn 0 0.0562
+## 13 location_Louisa -0.160 0.0562
+## 14 gender_male 1.07 0.0562
+## 15 frame_medium -0.320 0.0562
+## 16 frame_small -0.530 0.0562
+
+# plot variables ordered by importance (highest abs(coeff))
+model_final %>%
+extract_fit_parsnip() %>%
+vip(geom ="point") +
+theme_bw()
+
+
+
+
+Mean RMSE plus/minus standard error across repeated cross-validation folds as a function of lambda values
+
+
+
+
+
+
+
+Predicted BMI values against actual BMI values using final model for predicting test (unseen) data
+
+
+
+
+
+
+
+Top feature of importance, here measured as the features with highest absolute value of the Lasso regression coefficients from the final tuned model
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
\ No newline at end of file
diff --git a/session-tidymodels/case-study_files/figure-html/data-split-1.png b/session-tidymodels/docs/case-study_files/figure-html/data-split-1.png
similarity index 100%
rename from session-tidymodels/case-study_files/figure-html/data-split-1.png
rename to session-tidymodels/docs/case-study_files/figure-html/data-split-1.png
diff --git a/session-tidymodels/case-study_files/figure-html/load-data-1.png b/session-tidymodels/docs/case-study_files/figure-html/load-data-1.png
similarity index 100%
rename from session-tidymodels/case-study_files/figure-html/load-data-1.png
rename to session-tidymodels/docs/case-study_files/figure-html/load-data-1.png
diff --git a/session-tidymodels/case-study_files/figure-html/load-data-2.png b/session-tidymodels/docs/case-study_files/figure-html/load-data-2.png
similarity index 100%
rename from session-tidymodels/case-study_files/figure-html/load-data-2.png
rename to session-tidymodels/docs/case-study_files/figure-html/load-data-2.png
diff --git a/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-1.png b/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-1.png
new file mode 100644
index 0000000..3188bd1
Binary files /dev/null and b/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-1.png differ
diff --git a/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-2.png b/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-2.png
new file mode 100644
index 0000000..715bffe
Binary files /dev/null and b/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-2.png differ
diff --git a/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-3.png b/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-3.png
new file mode 100644
index 0000000..4371bca
Binary files /dev/null and b/session-tidymodels/docs/case-study_files/figure-html/unnamed-chunk-4-3.png differ
diff --git a/session-tidymodels/docs/index.html b/session-tidymodels/docs/index.html
index 6ff62b1..49330c9 100644
--- a/session-tidymodels/docs/index.html
+++ b/session-tidymodels/docs/index.html
@@ -108,7 +108,7 @@
