{"cells":[{"metadata":{},"cell_type":"markdown","source":"The objective of this notebook is to use simple linear models on the tabular data. I will try to keep it very simple to read and understand what is going on. Granted, these are not the best models, and I make no claim to have specific knowledge that would serve as a basis to these models. I just hope that the simplicity of the concepts in this notebook can give you some ideas that you could explore during these last days before the submissions deadline closes. But it can also be useful more generally by people that want to familiarise themselves with R.\n\nI have written before in the introduction of [another notebook](https://www.kaggle.com/douglaskgaraujo/bayesian-log-linear-model) why I am using log FVC values below."},{"metadata":{"_uuid":"051d70d956493feee0c6d64651c6a088724dca2a","_execution_state":"idle","trusted":true},"cell_type":"code","source":"library(tidyverse)\nlibrary(stargazer)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train <- read_csv(\"../input/osic-pulmonary-fibrosis-progression/train.csv\")\ntest <- read_csv(\"../input/osic-pulmonary-fibrosis-progression/test.csv\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"I created the function below to transform the dataset from the way we first load it into something closer I intend to use in the models. The reason I created the function instead of just applying these commands is that I want the exact same transformations to happen in the test data when we are preparing the submission file."},{"metadata":{"trusted":true},"cell_type":"code","source":"transformOriginalDataset <- function(dataset) {\n  # Uses log FVC to train the regression models\n  dataset$FVC_original <- dataset$FVC\n  dataset$FVC  <- log(dataset$FVC)\n  \n  # Transforms the 'Sex' variable into a dummy,\n  # a dummy for young patients, defined as those with less than 50 years of age,\n  # the smoking status variable into a yes/no currently smoker dummy,\n  # includes a variable with the patient-specific benchmark FVC\n  # (ie, FVC equal to 100% theoretical capacity),\n  # and finally includes the baseline Percent value\n  dataset <- dataset %>% \n    mutate(Male = ifelse(Sex == \"Male\", 1, 0),\n           CurrentlySmokes = ifelse(SmokingStatus == \"Currently smokes\", 1, 0),\n           YoungPatient = ifelse(Age < 50, 1, 0)) %>% \n    group_by(Patient) %>% \n    mutate(HundredPctFVC = 100 * mean(FVC_original / Percent),\n           FVCBase = head(FVC, 1),\n           PercentBase = head(Percent, 1))\n  \n  return(dataset)\n}","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"exclude_test_patient_data_from_trainset <- TRUE\n\nif(exclude_test_patient_data_from_trainset) {\n  train <- train %>% filter(!(Patient %in% unique(test$Patient)))\n}\n\ntrain <- rbind(train, test)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train <- transformOriginalDataset(train)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Below is the collection of models that we will be working with. As you can see, they increase in complexity from the 1st model to the 6th model, which has many more interaction terms. I will use model No 3 as an example below - mostly because its coefficients are essentially the same as simpler model No 1 but the interaction term allows me to illustrate to you how to include interaction terms in this notebook, in case you want to copy the notebook and  choose a different model to play around and submit a file correspondingly."},{"metadata":{"trusted":true},"cell_type":"code","source":"freq_pooled_models <- list()\nfreq_pooled_models[[1]] <- lm(FVC ~ Weeks + Male + CurrentlySmokes + YoungPatient + PercentBase, data = train)\nfreq_pooled_models[[2]] <- lm(FVC ~ Weeks + Male + CurrentlySmokes + YoungPatient + PercentBase + Weeks * Male, data = train)\nfreq_pooled_models[[3]] <- lm(FVC ~ Weeks + Male + CurrentlySmokes + YoungPatient + PercentBase + Weeks * CurrentlySmokes, data = train)\nfreq_pooled_models[[4]] <- lm(FVC ~ Weeks + Male + CurrentlySmokes + YoungPatient + PercentBase + Weeks * YoungPatient, data = train)\nfreq_pooled_models[[5]] <- lm(FVC ~ Weeks + Male + CurrentlySmokes + YoungPatient + PercentBase + Weeks * PercentBase, data = train)\nfreq_pooled_models[[6]] <- lm(FVC ~ Weeks + Male + CurrentlySmokes + YoungPatient + PercentBase + \n                                Weeks * CurrentlySmokes  +\n                                Weeks * Male  +\n                                Male * CurrentlySmokes  +\n                                Weeks * Male * CurrentlySmokes, data = train)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now we will use the excellent package ```stargazer```to compare the parameter estimaton of each of those models above."},{"metadata":{"trusted":true},"cell_type":"code","source":"stargazer(freq_pooled_models, type = \"text\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"How you pick one of the models, by selecting the index number (in this case, the 3rd model)."},{"metadata":{"trusted":true},"cell_type":"code","source":"model_used <- freq_pooled_models[[3]]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"If you customising this notebook, don't forget to confirm that you are including all coeficients for the model you chose. To do that, use the code below to check what the formula for that model is."},{"metadata":{"trusted":true},"cell_type":"code","source":"formula(model_used)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"So the output of ```formula(model_used)``` shows us that we should include: Weeks, Male, CurrentlySmokes, YoungPatient, PercentBase and the multiplication between Weeks and CurrentlySmokes. Let's see the corresponding coefficients (while also keeping tab of their \"names\" in the coefficient vector). The ```knitr::kable``` function is being used just to make it look better!\n\nOne thing to note here, which will be useful "},{"metadata":{"trusted":true},"cell_type":"code","source":"knitr::kable(coef(model_used))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"OK. So now that we have a model, let's use it to create predictions for each of those weeks that are required in the submission (ie, -12 weeks before the baseline CT and all the way until 133 weeks after). I tried making my code as verbose as possible so that you can understand what I am doing. Functions in the ```tidyverse``` have excellent documentation and loads of community questions & answers, so don't hesitate to investigate if you have not understood any specific point."},{"metadata":{"trusted":true},"cell_type":"code","source":"test_sub <- transformOriginalDataset(test)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Variable ```test_sub``` has now undergone the transformations that were also applied to the training test when used it to model FVC outcomes. But it is not yet ready to be used to construct the submission file. To understand why, note that the submission rules require us to generate a predicted FVC and confidence for each patient in the test set during a large range of weeks (from -12 to 133, as I'm sure you know). So now we need some supporting variables and a code that will create such a data.frame for us, and then fill it out with the predictions for FVC and confidence values."},{"metadata":{"trusted":true},"cell_type":"code","source":"weeks <- seq(from = -12, to = 133, by = 1)\nPatient_Week <- paste(rep(unique(test$Patient), each = length(weeks)), \n                      seq(from = -12, to = 133, by = 1), sep = \"_\")\n\ntest_sub <- test_sub %>% \n  group_by(Patient) %>% \n  summarise(FVCBase = head(FVCBase, 1),\n            Male = head(Male, 1),\n            CurrentlySmokes = head(CurrentlySmokes, 1),\n            YoungPatient = head(YoungPatient, 1),\n            PercentBase = head(PercentBase, 1)) %>% \n  uncount(length(weeks)) %>% \n  cbind(Patient_Week, Weeks = rep(seq(from = -12, to = 133, by = 1), length(unique(test$Patient))), .)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"In the first stage above, we merely compiled the patient-specific information we will use to predict the values. Note that each patient only has one value here. The reason is that our predictions will be measured against test dataset containing only one row of information per patient. So our prediction needs to be based on the first value only as well. That is what we will be doing next. Just remember to keep it consistent with ```formula(model_used)``` if you adapt this notebook."},{"metadata":{"trusted":true},"cell_type":"code","source":"test_sub <- test_sub %>% \n  mutate(FVC = exp(coef(model_used)[which(names(coef(model_used)) == \"(Intercept)\")] +\n                   coef(model_used)[which(names(coef(model_used)) == \"Weeks\")] * Weeks +\n                   coef(model_used)[which(names(coef(model_used)) == \"Male\")] * Male +\n                   coef(model_used)[which(names(coef(model_used)) == \"CurrentlySmokes\")] * CurrentlySmokes +\n                   coef(model_used)[which(names(coef(model_used)) == \"YoungPatient\")] * YoungPatient +\n                   coef(model_used)[which(names(coef(model_used)) == \"Weeks:CurrentlySmokes\")] * Weeks * CurrentlySmokes\n                   ),\n         Confidence = (exp(1.96 * summary(model_used)$sigma) - 1) * FVC) %>% \n  select(Patient_Week, FVC, Confidence)\n\nhead(test_sub)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Great! Now let's write the submission file."},{"metadata":{"trusted":true},"cell_type":"code","source":"write_csv(test_sub, \"submission.csv\")","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"name":"ir","display_name":"R","language":"R"},"language_info":{"name":"R","codemirror_mode":"r","pygments_lexer":"r","mimetype":"text/x-r-source","file_extension":".r","version":"3.6.3"}},"nbformat":4,"nbformat_minor":4}