---
title: "Uppstartslektion Modul 4"
format: html
editor: visual
editor_options: 
  chunk_output_type: console
---

Nedan finns ett förslag på arbetsgång för att bekanta sig med datamaterialet och för att påbörja letandet efter den "bästa" modellen. Vad som är den bästa modellen beror på vilka kriterier som används. Vi har alltså inget facit. Er uppgift är att övertyga oss om varför er modell är den bästa. Alla grupper har olika datamaterial vilket gör att ni kommer att få olika modeller.

```{r}
#| message: false
# Ladda R-paket
library(tidyverse)
library(GGally)
library(gtsummary)
library(ggrepel)
library(rsample)
```

1.  Öppna filen med datamaterialet ni ska arbeta med. Filen finns att hämta på sidan *Modul 4: Om modul 4 (Instruktioner, datamaterial m.m.)* i Canvas. Välj en responsvariabel och **ta bort de tre responsvariabler som inte ska användas** så ni inte av misstag använde de som förklaringsvariabler. **Filtera ut den karosstyp er grupp blivit tilldelad** (se dokumentet "Data uppdelat på grupp inkl kod till hjälp" på Canvas).

    ```{r}
    # Öppna filen
    bildata <- readRDS("bildata.rds")
    ```

    ```{r}
    # Ta bort de responsvariabler som inte ska användas

    # Ex. om ni vill ha koldioxidutsläpp som respons:
    bildata <- 
      bildata %>% 
      select(-c(Bransle_blandad, Bransle_stad, Fordonsskatt, Bransle_landsvag))
    ```

    ```{r}
    # Vad heter de olika karosstyperna?
    unique(bildata$BODY)

    # Ex. om ni ska ha karosstyp SEDAN:
    bildata <- 
      bildata %>% 
      filter(BODY == "KOMBI-SEDAN") # Se canvassida för information om er grupps karosstyp

    ```

2.  **Dela upp datamaterialet i träningsdata och testdata med hjälp av ert tilldelade frö (seed)**. Se dokumentet "Data uppdelat på grupp inkl kod till hjälp" på Canvas. Spara sedan datamaterialen.

    ```{r}
    # Dela upp datamaterialet (80% till träning och 20% till utvärdering)

    set.seed(687) # Se canvassida för information om er grupps slumptalsfrö

    bildata_split <- initial_split(bildata, prop = 0.80)
    bildata_train <- training(bildata_split)
    bildata_test <- testing(bildata_split)
    ```

    ```{r}
    # Spara datamaterialen

    saveRDS(bildata_train, file = "bildata_train_grupp_x.rds")
    saveRDS(bildata_test, file = "bildata_test_grupp_x.rds")
    ```

3.  Lär känna varje enskild variabel genom att t.ex. göra lådagram för de kvantitativa variablerna. Har någon kvantitativ variabel outliers? Vilka variabler är kvalitativa? Innehåller någon variabel saknade värden? Använd gärna verktyget **ggplot2 builder** i Addins. *Se Lektion 1 för instruktioner.*

    ```{r}
    ## Klistra in kod från ggplot2 builder här

    ggplot(bildata_test) +
       aes(x = Koldioxidutslapp, y = Effekt) +
       geom_point(colour = "#112446") +
       geom_smooth(
           span = 1L,
           level = 1L,
           method = "lm",
           se = FALSE,
           colour = "#112446"
       ) +
       theme_minimal()
    ```

4.  Studera sambandet mellan par av kvantitativa variabler genom att plotta numeriska variabler parvis mot varandra och beräkna korrelationer. *Se Lektion 7.*

    ```{r}
    bildata_train %>% 
      select(where(is.numeric)) %>% # Välj ut numeriska variabler
      ggpairs()
    ```

5.  Rita ett spridningsdiagram mellan responsvariabeln och den förklaringsvariabel som ni tror har störst betydelse för att förklara er responsvariabel. Hur ser sambandet ut? *Använd **ggplot 2 builder** eller skriv egen kod*. *Se Lektion 1 för instruktioner.*

    ```{r}
    # Lägg in ggplot-kod nedan
    lm()
    ```

6.  Gör en enkel linjär regressionsanalys med den förklaringsvariabel ni valde i uppgift 5. Vilken förklaringsgrad har modellen? Vilken justerad förklaringsgrad har modellen? *Se Lektion 7*.

    ```{r}
    Modell_1 <- 
      bildata_train %>% 
      lm(responsvariabel ~ förklaringsvariabel,  # ändra variabelnamn
         data = .) 
      
    Modell_1 %>% summary()
      
    ```

7.  Gör om spridningsdiagrammet från uppgift 5 med tillägget att punkternas färg beror på bränsletyp (bensin/diesel). Ger diagrammet indikation på att den kvalitativa variabeln bränsletyp bör tas med i regressionsanalysen? *Se Lektion 8*.

    ```{r}
    # Klistra in ggplot-kod nedan
    bildata_train %>%
      ggplot(
        aes(
          x = förklaringsvariabel, # ändra variabelnamn
          y = responsvariabel, # ändra variabelnamn
          color = Bransletyp
        )
      ) +
      geom_point() +
      geom_smooth(
        method = "lm",
        se = FALSE
      )
      
    ```

8.  Gör en ny regressionsanalys som också inkluderar bränsletyp som förklarande variabel. Har bränsletyp ett signifikant samband (5% signifikansnivå) med responsvariabeln, givet att vi att vi har med den andra förklaringsvariabeln i modellen? *Se Lektion 8*.

    ```{r}
    Modell_2 <- 
      bildata_train %>% 
      lm(responsvariabel ~ förklaringsvariabel + Bransletyp, # ändra variabelnamn
         data = .) 

    Modell_2 %>% summary()
    ```

9.  Hur förändras $R^2$ och $R^2_{adj}$ jämfört med modellen i uppgift 6? När används med fördel den justerade förklaringsgraden framför den "vanliga" förklaringsgraden, $R^2$?

10. Vilken modell har lägst prediktions-MSE? Stämmer detta överens med din slutsats från uppgift 9 baserat på $R^2_{adj}$ värdena? *Se Lektion 7.*

    ```{r}

    #Prediktera på testdata
    pred1 <- predict(Modell_1, newdata = bildata_test)
    pred2 <- predict(Modell_2, newdata = bildata_test)

    #Beräkna MSE
    MSE1 = mean((bildata_test$responsvariabel - pred1)^2,
                na.rm = TRUE) # ändra variabelnamn
    MSE2 = mean((bildata_test$responsvariabel - pred2)^2,
                na.rm = TRUE) # ändra variabelnamn

    print(paste("Prediktions MSE Modell 1:", round(MSE1)))
    print(paste("Prediktions MSE Modell 2:", round(MSE2)))
    ```

11. Gör en residualanalys för modellen i uppgift 8. *Se Lektion 7*. Vilka bilar har störst residualer? Finns det någon gemensam nämnare med dessa bilar?

    ```{r}
    # Residualplottar: 

    # Residuals vs fitted 
      plot(Modell_2, which=1)
    # QQ-plot
      plot(Modell_2, which=2)
    ```

    ```{r}
    # Här är ett alternativ om du vill göra snyggare residauls vs fitted-plot med ggplot genom att först spara "predikterade värden" och "residualer" som nya variabler  i ert datamaterial och sedan använda ggplot

    # Lägg till predikterade värden i bildata_train
    bildata_train$y_hat <- predict(Modell_2, newdata = bildata_train)

    # Lägg till residualer i bildata_train
    bildata_train <- bildata_train %>% mutate(res = responsvariabel - y_hat)  # ändra variabelnamn

    # Residuals vs fitted-plot med ggplot
    bildata_train %>%
      ggplot(aes(x = y_hat,
                 y = res
                 )) +
      geom_point() +
      labs(
        x = "Fitted values",
        y = "Residuals") +
      theme_minimal()

    # Residuals vs fitted-plot med ggplot där text med bilmärke och version anges för de residualser som är större än 80 eller mindre än -80 (detta värde kan ändras)
    bildata_train %>% 
      ggplot(aes(x = y_hat,
                 y = res
                 )) +
      geom_point() +
      labs(
        x = "Fitted values",
        y = "Residuals") +
      geom_text_repel(
        data = filter(bildata_train, abs(res) > 80), # ändra värde för fler/färre
        aes(label = paste(Bilmarke, Version))
                      ) +
      theme_minimal()

    ```

12. Finns det indikationer på att en samspelsterm mellan bränsletyp och den kvantitativa variabeln som ni har i modellen i uppgift 8 bör inkluderas? Ett sätt att ta reda på det är att studera spridningsdiagrammet i uppgift 7. Är de skattade linjerna parallella? Om det finns indikationer på att ni behöver en samspelsterm, testa då om samspelet är signifikant (dvs. om effekten som den kvantitativa variabeln har på responsen skiljer sig mellan bensin- och dieselbilar). *Se Lektion 8*.

    ```{r}
    Modell_3 <- 
      bildata_train %>% 
      lm(responsvariabel ~ förklaringsvariabel * Bransletyp,  # ändra variabelnamn
          
         data = .) 

    Modell_3 %>% summary()
    ```

13. Prova att utöka modellen med ytterligare en variabel. Blev modellen bättre? Vilka mått är intressanta att titta på då modeller ska jämföras?

14. Fortsätt att undersöka om fler variabler bör tas med i modellen.

**Notera**:\
Att bara utöka modellen med fler och fler variabler leder inte alltid till den ”bästa” modellen. Den variabel som är starkast korrelerad med responsvariabeln är den variabel som är “bäst” om vi bara får välja en variabel i modellen. I modellbyggnadsprocessen gäller det att hitta den kombination av variabler som tillsammans förklarar så mycket som möjligt av variationen i responsvariabeln, dvs. den kombination av variabler som kompletterar varandra så bra som möjligt. Den kombinationen behöver inte nödvändigtvis innehålla den variabeln som är starkast korrelerad med responsvariabeln.

Notera även att det inte nödvändigtvis är den modell med högst $R^2_{adj}$ som är bäst på prediktion. Det kan hända att ni har överanpassat modellen till ert datamaterial vilket då gör den sämre för prediktion av andra bilar (testdatat).

Att komma fram till den bästa modellen är inte lätt. Det är en hel del testande fram och tillbaka. När ni har kommit fram till den modell som ni anser är den bästa ska ni börja skriva rapporten. Fundera över vilken strategi/vilket tillvägagångssätt ni har använt för att komma fram till den modell som ni anser är den bästa. Den strategin ska beskrivas i er rapport.

**Läs rapportinstruktionerna och titta på modelltexten innan ni börjar skriva!**
