Raport 3-Regresja nieparametryczna

library(car)
library(KernSmooth)
library(tidyverse)
library(MASS)
data("Prestige")
attach(Prestige)
data("mcycle")

Prestige data

Zbiór danych przedstawia wskaźnik presitiżu zawodu według skali Pineo-Portera na podstawie badania społecznego przeprowadzonego w połowie lat 60. W pozostałych kolumnach znajdują się informacje dotyczące średniego wykszatłcenia, dochodu, tego jaki procent pracujących w tym zawodzie stanowią kobiety oraz przyporządkowania danego zawodu do jednego z trzech typów: Professional/Managerial, Blue Collar lub White Collar.

Przykład 1: Prestiż względem dochodu. Zestaw danych “Prestige” (z pakietu “car”) zawiera dane nt. prestiżu n=102 Kanadyjskich zawodów z 1971 roku, a także średni dochód w danym zawodzie. Do zbadania zależności między prestiżem a dochodem wykorzystaj metody regresji nieparametrycznej.

# Relacja między dochodem a prestiżem
ggplot(Prestige) +
  geom_point(aes(x=income,y=prestige))

local polynomial regression- metoda lokalnego wygładzania

Główne parametry, takie jak degree i bandwidth, mają kluczowy wpływ na sposób, w jaki wygładzana jest zależność.

Parametr degree określa stopień wielomianu użytego do lokalnej aproksymacji. Degree=0 oznacza użycie lokalnych średnich (co odpowiada regresji najbliższych sąsiadów). Degree=1 oznacza wykorzystanie lokalnej regresji liniowej Degree=2 oznacza wykorzystanie lokalnej regresji kwadratowej.

Parametr bandwidth decyduje o rozpiętości jądra, czyli definiuje rozmiar obszaru w jakim obserwacje wpływają na funkcję gęstości.Im większa wartość bandwidth, tym bardziej wygładzona jest funkcja. Zarówno zbyt niska wartość bandwidth może być nieodpowiednia (overfitting- model zbyt szczegółowy), jak i również zbyt wysoka wartość tego parametru (underfitting- zbyt duże wygładzenie, zbyt duże uogólnienie).

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
# lokalna regresja wielomianowa, gdzie dregree=0 i bandwitch=5
fit <- locpoly(prestige,income,
        degree=0, bandwidth=5) %>% as_tibble

ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=fit, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
# lokalna regresja wielomianowa, gdzie dregree=1 i bandwitch=5
fit2 <- locpoly(prestige, income,
        degree=1, bandwidth=5) %>% as_tibble

ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=fit2, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
# lokalna regresja wielomianowa, gdzie dregree=2 i bandwitch=5
fit3 <- locpoly(prestige, income,
        degree=2, bandwidth=5) %>% as_tibble

ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=fit3, aes(x=x,y=y), col='blue')

#Degree Można zauważyć, że wybór degree wpływa na funkcję gęstości na wykresie. W przypadku degree=2 występuje overfitting. Lepszym wyborem byłoby degree równe 0 lub 1.

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
# lokalna regresja wielomianowa, gdzie dregree=0 i bandwitch=2
fit4 <- locpoly(prestige, income,
        degree=0, bandwidth=2) %>% as_tibble

ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=fit4, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
# lokalna regresja wielomianowa, gdzie dregree=0 i bandwitch=5
fit <- locpoly(prestige, income,
        degree=0, bandwidth=5) %>% as_tibble

ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=fit, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
# lokalna regresja wielomianowa, gdzie dregree=0 i bandwitch=10
fit5 <- locpoly(prestige, income,
        degree=0, bandwidth=10) %>% as_tibble

ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=fit5, aes(x=x,y=y), col='blue')

#Bandwidth W przypadku bandwidth=2 funkcja jest zbyt dopasowana (overfitting), zaś przy bandwidth=10 powstaje zbyt uogólniona postać funkcji. Najlepszym wyborem byłoby bandwidth=5.

Smooth.spline- Sploty interpolujące

Parametr CV (cross-validation) pozwala na automatyczny wybór pasma.

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
fit6 <-smooth.spline(prestige, income, cv=TRUE) 

smr <- data.frame(x=fit6$x,y=fit6$y)
ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=smr, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
fit7 <-smooth.spline(prestige, income, spar=0.2) 

smr2 <- data.frame(x=fit7$x,y=fit7$y)
ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=smr2, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
fit8 <-smooth.spline(prestige, income, spar=0.9) 

smr3 <- data.frame(x=fit8$x,y=fit8$y)
ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=smr3, aes(x=x,y=y), col='blue')

# zależność pomiędzy prestiżem zawodu a wysokością dochodu
fit9 <-smooth.spline(prestige, income, spar=0.75) 

smr4 <- data.frame(x=fit9$x,y=fit9$y)
ggplot(Prestige) +
  geom_point(aes(x=prestige,y=income)) +
  geom_line(data=smr4, aes(x=x,y=y), col='blue')

Spar=0.2 powoduje zbyt mocne dopasowanie do danych. Spar=0.9 prezentuje zbyt duże uogólnienie. Wartość parametru w około 0.75 wydaje się odpowiednia, wykres prezentuje się podobnie do wybranego automatycznie dzięki parametrowi CV (cross-validation).

mcycle data

Przykład 2: Wypadek motocyklowy. Zbiór danych “mcycle” (z pakietu MASS) zawiera n=133 pary punktów czasowych (w ms) i obserwowanych przyspieszeń głowy (w g), które zostały zarejestrowane w symulowanym wypadku motocyklowym. Do zbadania zależności między czasem a przyspieszeniem wykorzystaj metody regresji nieparametrycznej.

# zależność między czasem a przyspieszeniem
ggplot(mcycle) +
  geom_point(aes(x=times,y=accel))

local polynomial regression- metoda lokalnego wygładzania

# zależność pomiędzy czasem a przyspieszeniem
# lokalna regresja wielomianowa, gdzie dregree=0 i bandwitch=5
fit10 <- locpoly(mcycle$times, mcycle$accel,
        degree=0, bandwidth=5) %>% as_tibble

ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=fit10, aes(x=x,y=y), col='blue')

# zależność pomiędzy czasem a przyspieszeniem
# lokalna regresja wielomianowa, gdzie dregree=1 i bandwitch=5
fit11 <- locpoly(mcycle$times, mcycle$accel,
        degree=1, bandwidth=5) %>% as_tibble

ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=fit11, aes(x=x,y=y), col='blue')

# zależność pomiędzy czasem a przyspieszeniem
# lokalna regresja wielomianowa, gdzie dregree=2 i bandwitch=5
fit12 <- locpoly(mcycle$times, mcycle$accel,
        degree=2, bandwidth=5) %>% as_tibble

ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=fit12, aes(x=x,y=y), col='blue')

#Degree W przypadku zbioru danych dotyczącego wypadków motocyklowych wybór dregree=0 daje zbyt ogólną postać fukncji. Najbardziej dopasowaną do tych danych funkcją wydaje się być ta wykorzystująca lokalną funkcję kwadratową (degree=2).

# zależność pomiędzy czasem a przyspieszeniem
# lokalna regresja wielomianowa, gdzie dregree=2 i bandwitch=2
fit13 <- locpoly(mcycle$times, mcycle$accel,
        degree=2, bandwidth=2) %>% as_tibble

ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=fit13, aes(x=x,y=y), col='blue')

# zależność pomiędzy czasem a przyspieszeniem
# lokalna regresja wielomianowa, gdzie dregree=2 i bandwitch=10
fit14 <- locpoly(mcycle$times, mcycle$accel,
        degree=2, bandwidth=10) %>% as_tibble

ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=fit14, aes(x=x,y=y), col='blue')

#Bandwidth Szerokość pasma 5 jest najbardziej odpowiednia, gdyż bandwidth=2 zbyt szczegółowo przypasowuje funkcję do danych, zaś bandwidth=10 tworzy zbyt ogólną funkcję.

Smooth.spline- Sploty interpolujące

# zależność pomiędzy czasem a przyspieszeniem
fit15 <-smooth.spline(mcycle$times, mcycle$accel, cv=TRUE) 

smr5 <- data.frame(x=fit15$x,y=fit15$y)
ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=smr5, aes(x=x,y=y), col='blue')

# zależność pomiędzy czasem a przyspieszeniem
fit16 <-smooth.spline(mcycle$times, mcycle$accel, spar=0.65) 

smr6 <- data.frame(x=fit16$x,y=fit16$y)
ggplot(mcycle) +
  geom_point(aes(x=times,y=accel)) +
  geom_line(data=smr6, aes(x=x,y=y), col='blue')

Przy spar= 0.65 funkcja jest odpowiednio wygładzona. Niższa wartość parametru mogłaby zbyt szczegółówo dopasować krzywą do danych.