Regresja i metody uczenia nadzorowanego - Prognozowanie cen nieruchomości
Jakub Miros, Grzegorz Pakier
2026-06-19
Wprowadzenie
Celem niniejszego projektu jest budowa zaawansowanych modeli uczenia nadzorowanego dedykowanych problemowi regresji, jakim jest wycena wartości nieruchomości. W ramach analizy zweryfikowana zostanie skuteczność klasycznych podejść regresyjnych w zderzeniu z modelami nieliniowymi i zespołowymi. Końcowym efektem prac jest wybór optymalnego modelu, który charakteryzuje się najwyższą zdolnością generalizacji i najniższym poziomem błędu predykcji na danych testowych.
1. Eksploracja danych
1.1. Inicjacja i wczytanie zbioru danych
W pierwszym kroku ładowane są niezbędne pakiety narzędziowe odpowiedzialne za przetwarzanie danych, wizualizację, modelowanie oraz ekstrakcję metryk jakościowych.
library(tidyverse)
library(caret)
library(corrplot)
library(skimr)
library(DataExplorer)
library(rpart)
library(rpart.plot)
library(partykit)
library(randomForest)
library(ranger)
library(neuralnet)
library(NeuralNetTools)
# Wylaczenie notacji naukowej i ustawienie lokalizacji komunikatow
options(scipen = 999)
Sys.setenv(LANG = "pl")
# Wczytanie danych surowych
df_houses <- read.csv("r5_houses.csv", stringsAsFactors = FALSE)1.2 Struktura i statystyki opisowe
Przeprowadzenie wstępnego wglądu w strukturę ramki danych pozwala zidentyfikować typy zmiennych (numeryczne ciągłe, dyskretne oraz potencjalne zmienne kategoryczne zapisane jako liczby) oraz poznać podstawowe charakterystyki rozkładów.
cat("--- Struktura zbioru danych ---\n")## --- Struktura zbioru danych ---glimpse(df_houses)## Rows: 21,408
## Columns: 22
## $ price <dbl> 221900, 538000, 180000, 604000, 510000, 1225000, 257500,…
## $ bedrooms <int> 3, 3, 2, 4, 3, 4, 3, 3, 3, 3, 3, 2, 3, 3, 5, 4, 3, 4, 2,…
## $ bathrooms <dbl> 1.00, 2.25, 1.00, 3.00, 2.00, 4.50, 2.25, 1.50, 1.00, 2.…
## $ sqft_living <int> 1180, 2570, 770, 1960, 1680, 5420, 1715, 1060, 1780, 189…
## $ sqft_lot <int> 5650, 7242, 10000, 5000, 8080, 101930, 6819, 9711, 7470,…
## $ floors <dbl> 1.0, 2.0, 1.0, 1.0, 1.0, 1.0, 2.0, 1.0, 1.0, 2.0, 1.0, 1…
## $ waterfront <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ view <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0,…
## $ condition <int> 3, 3, 3, 5, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 4, 4,…
## $ grade <int> 7, 7, 6, 7, 8, 11, 7, 7, 7, 7, 8, 7, 7, 7, 7, 9, 7, 7, 7…
## $ sqft_above <int> 1180, 2170, 770, 1050, 1680, 3890, 1715, 1060, 1050, 189…
## $ sqft_basement <int> 0, 400, 0, 910, 0, 1530, 0, 0, 730, 0, 1700, 300, 0, 0, …
## $ yr_built <int> 1955, 1951, 1933, 1965, 1987, 2001, 1995, 1963, 1960, 20…
## $ yr_renovated <int> 0, 1991, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ zipcode <int> 98178, 98125, 98028, 98136, 98074, 98053, 98003, 98198, …
## $ lat <dbl> 47.5112, 47.7210, 47.7379, 47.5208, 47.6168, 47.6561, 47…
## $ long <dbl> -122.257, -122.319, -122.233, -122.393, -122.045, -122.0…
## $ house_age <int> 60, 64, 82, 50, 28, 14, 20, 52, 55, 12, 50, 73, 88, 38, …
## $ was_renovated <int> 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ total_sqft <int> 1180, 2970, 770, 2870, 1680, 6950, 1715, 1060, 2510, 189…
## $ living_ratio <dbl> 0.20881260, 0.35482535, 0.07699230, 0.39192162, 0.207895…
## $ is_luxury <int> 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…cat("\n--- Podstawowe statystyki opisowe ---\n")##
## --- Podstawowe statystyki opisowe ---summary(df_houses)## price bedrooms bathrooms sqft_living
## Min. : 75000 Min. : 0.000 Min. :0.000 Min. : 290
## 1st Qu.: 320000 1st Qu.: 3.000 1st Qu.:1.500 1st Qu.:1420
## Median : 450000 Median : 3.000 Median :2.250 Median :1900
## Mean : 518937 Mean : 3.362 Mean :2.098 Mean :2052
## 3rd Qu.: 635000 3rd Qu.: 4.000 3rd Qu.:2.500 3rd Qu.:2520
## Max. :1999950 Max. :33.000 Max. :7.500 Max. :7730
## sqft_lot floors waterfront view
## Min. : 520 Min. :1.000 Min. :0.000000 Min. :0.0000
## 1st Qu.: 5027 1st Qu.:1.000 1st Qu.:0.000000 1st Qu.:0.0000
## Median : 7578 Median :1.500 Median :0.000000 Median :0.0000
## Mean : 15036 Mean :1.491 Mean :0.005325 Mean :0.2168
## 3rd Qu.: 10556 3rd Qu.:2.000 3rd Qu.:0.000000 3rd Qu.:0.0000
## Max. :1651359 Max. :3.500 Max. :1.000000 Max. :4.0000
## condition grade sqft_above sqft_basement
## Min. :1.000 Min. : 1.000 Min. : 290 Min. : 0.0
## 1st Qu.:3.000 1st Qu.: 7.000 1st Qu.:1190 1st Qu.: 0.0
## Median :3.000 Median : 7.000 Median :1550 Median : 0.0
## Mean :3.408 Mean : 7.628 Mean :1767 Mean : 284.6
## 3rd Qu.:4.000 3rd Qu.: 8.000 3rd Qu.:2190 3rd Qu.: 550.0
## Max. :5.000 Max. :13.000 Max. :7420 Max. :3260.0
## yr_built yr_renovated zipcode lat
## Min. :1900 Min. : 0.00 Min. :98001 Min. :47.16
## 1st Qu.:1951 1st Qu.: 0.00 1st Qu.:98033 1st Qu.:47.47
## Median :1975 Median : 0.00 Median :98065 Median :47.57
## Mean :1971 Mean : 82.23 Mean :98078 Mean :47.56
## 3rd Qu.:1997 3rd Qu.: 0.00 3rd Qu.:98118 3rd Qu.:47.68
## Max. :2015 Max. :2015.00 Max. :98199 Max. :47.78
## long house_age was_renovated total_sqft
## Min. :-122.5 Min. : 0.00 Min. :0.0000 Min. : 290
## 1st Qu.:-122.3 1st Qu.: 18.00 1st Qu.:0.0000 1st Qu.:1490
## Median :-122.2 Median : 40.00 Median :0.0000 Median :2160
## Mean :-122.2 Mean : 44.01 Mean :0.0412 Mean :2336
## 3rd Qu.:-122.1 3rd Qu.: 64.00 3rd Qu.:0.0000 3rd Qu.:2960
## Max. :-121.3 Max. :115.00 Max. :1.0000 Max. :9950
## living_ratio is_luxury
## Min. :0.00061 Min. :0.00000
## 1st Qu.:0.15622 1st Qu.:0.00000
## Median :0.24701 Median :0.00000
## Mean :0.32329 Mean :0.06806
## 3rd Qu.:0.40658 3rd Qu.:0.00000
## Max. :4.64491 Max. :1.00000# Sprawdzenie ostatecznych wymiarow zbioru danych
cat("\nWymiary zbioru (wiersze, kolumny):\n")##
## Wymiary zbioru (wiersze, kolumny):dim(df_houses)## [1] 21408 221.3 Analiza braków danych i wariancji predyktorów
Przed przystąpieniem do modelowania konieczne jest zweryfikowanie kompletności bazy danych pod kątem występowania brakujących wartości oraz eliminacja zmiennych bezwartościowych, czyli takich, których wariancja jest bliska zeru.
cat("--- Liczba braków danych w kolumnach ---\n")## --- Liczba braków danych w kolumnach ---missing_values <- colSums(is.na(df_houses))
print(missing_values)## price bedrooms bathrooms sqft_living sqft_lot
## 0 0 0 0 0
## floors waterfront view condition grade
## 0 0 0 0 0
## sqft_above sqft_basement yr_built yr_renovated zipcode
## 0 0 0 0 0
## lat long house_age was_renovated total_sqft
## 0 0 0 0 0
## living_ratio is_luxury
## 0 0cat("\n--- Analiza zmiennych near-zero variance ---\n")##
## --- Analiza zmiennych near-zero variance ---nzv_analysis <- nearZeroVar(df_houses, saveMetrics = TRUE)
print(nzv_analysis)## freqRatio percentUnique zeroVar nzv
## price 1.000000 18.226831091 FALSE FALSE
## bedrooms 1.441826 0.060724963 FALSE FALSE
## bathrooms 1.392264 0.126121076 FALSE FALSE
## sqft_living 1.022222 4.591741405 FALSE FALSE
## sqft_lot 1.238754 45.184043348 FALSE FALSE
## floors 1.312778 0.028026906 FALSE FALSE
## waterfront 186.789474 0.009342302 FALSE TRUE
## view 20.770053 0.023355755 FALSE TRUE
## condition 2.468750 0.023355755 FALSE FALSE
## grade 1.481115 0.056053812 FALSE FALSE
## sqft_above 1.009524 4.264760837 FALSE FALSE
## sqft_basement 59.684932 1.363976084 FALSE TRUE
## yr_built 1.239910 0.541853513 FALSE FALSE
## yr_renovated 225.560440 0.326980568 FALSE TRUE
## zipcode 1.020339 0.326980568 FALSE FALSE
## lat 1.000000 23.514573991 FALSE FALSE
## long 1.027027 3.512705531 FALSE FALSE
## house_age 1.239910 0.541853513 FALSE FALSE
## was_renovated 23.272109 0.009342302 FALSE TRUE
## total_sqft 1.008065 5.222346786 FALSE FALSE
## living_ratio 1.200000 94.282511211 FALSE FALSE
## is_luxury 13.693205 0.009342302 FALSE FALSE1.4 Analiza korelacji zmiennych numerycznych
W celu zbadania współliniowości cech ciągłych wyznaczono macierz korelacji liniowej Pearsona. Identyfikacja silnie skorelowanych predyktorów pozwala na uniknięcie problemów numerycznych i błędów interpretacyjnych w modelach strukturalnych.
# Do korelacji wybieramy tylko zmienne numeryczne
numeric_cols <- df_houses %>%
select(where(is.numeric))
# Obliczenie macierzy korelacji
data_corr <- cor(numeric_cols, use = "complete.obs")
# Wyświetlenie macierzy korelacji
corrplot(data_corr,
method = "color",
type = "lower",
order = "hclust",
tl.col = "black",
tl.srt = 45,
addCoef.col = "black",
number.cex = 0.6,
diag = TRUE,
title = "Korelacje zmiennych numerycznych",
mar = c(0,0,1,0))
2. Przygotowanie danych
2.1 Czyszczenie zmiennych i transformacja typów
Zmienne jakościowe i dychotomiczne zostały przekształcone na typ faktorowy. Ze względu na fakt, że najwyższa klasa jakości wykonania budynku (grade = 13) reprezentuje skrajnie małą liczbę obserwacji i mogłaby nie wylosować się do zbioru treningowego (generując błędy braku dopasowania poziomów w modelach regresji liniowej oraz sieciach neuronowych - co było problemem w przy pierwszej próbie), została ona na etapie czyszczenia połączona z klasą 12. Następnie zweryfikowano rozkłady cech za pomocą histogramów.
df_cleaned <- df_houses %>%
mutate(
# Polaczenie rzadkiej klasy grade 12 i 13 w celu stabilizacji estymacji
grade = ifelse(grade == 13, 12, grade),
# Konwersja zmiennych kategorycznych na czynniki
waterfront = as.factor(waterfront),
view = as.factor(view),
condition = as.factor(condition),
grade = as.factor(grade),
was_renovated = as.factor(was_renovated),
is_luxury = as.factor(is_luxury),
zipcode = as.factor(zipcode)
)
cat("--- Końcowa struktura danych po modyfikacji typów ---\n")## --- Końcowa struktura danych po modyfikacji typów ---str(df_cleaned)## 'data.frame': 21408 obs. of 22 variables:
## $ price : num 221900 538000 180000 604000 510000 ...
## $ bedrooms : int 3 3 2 4 3 4 3 3 3 3 ...
## $ bathrooms : num 1 2.25 1 3 2 4.5 2.25 1.5 1 2.5 ...
## $ sqft_living : int 1180 2570 770 1960 1680 5420 1715 1060 1780 1890 ...
## $ sqft_lot : int 5650 7242 10000 5000 8080 101930 6819 9711 7470 6560 ...
## $ floors : num 1 2 1 1 1 1 2 1 1 2 ...
## $ waterfront : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ view : Factor w/ 5 levels "0","1","2","3",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ condition : Factor w/ 5 levels "1","2","3","4",..: 3 3 3 5 3 3 3 3 3 3 ...
## $ grade : Factor w/ 11 levels "1","3","4","5",..: 6 6 5 6 7 10 6 6 6 6 ...
## $ sqft_above : int 1180 2170 770 1050 1680 3890 1715 1060 1050 1890 ...
## $ sqft_basement: int 0 400 0 910 0 1530 0 0 730 0 ...
## $ yr_built : int 1955 1951 1933 1965 1987 2001 1995 1963 1960 2003 ...
## $ yr_renovated : int 0 1991 0 0 0 0 0 0 0 0 ...
## $ zipcode : Factor w/ 70 levels "98001","98002",..: 67 56 17 59 38 30 3 69 61 24 ...
## $ lat : num 47.5 47.7 47.7 47.5 47.6 ...
## $ long : num -122 -122 -122 -122 -122 ...
## $ house_age : int 60 64 82 50 28 14 20 52 55 12 ...
## $ was_renovated: Factor w/ 2 levels "0","1": 1 2 1 1 1 1 1 1 1 1 ...
## $ total_sqft : int 1180 2970 770 2870 1680 6950 1715 1060 2510 1890 ...
## $ living_ratio : num 0.209 0.355 0.077 0.392 0.208 ...
## $ is_luxury : Factor w/ 2 levels "0","1": 1 1 1 1 1 2 1 1 1 1 ...# Generowanie histogramów dla zbioru surowego
plot_histogram(
data = df_houses %>% select(where(is.numeric)),
title = "Histogramy zmiennych numerycznych (Zbiór surowy)",
ggtheme = theme_minimal(),
theme_config = list(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14)
)
)
2.2 Podział danych i weryfikacja stratyfikacji
Zbiór danych został podzielony na część treningową (70%) i testową (30%) z zachowaniem stratyfikacji względem zmiennej celu (price). Identyczność rozkładów w obu podzbiorach zweryfikowano za pomocą statystyk opisowych oraz wykresu gęstości.
set.seed(123)
train_indices <- createDataPartition(df_cleaned$price, p = 0.7, list = FALSE)
housing.train <- df_cleaned[train_indices, ]
housing.test <- df_cleaned[-train_indices, ]
# Wyciągnięcie wektorów cen w celu kontroli rozkładów
real.train <- housing.train$price
real.test <- housing.test$price
# Definicja funkcji do wyciągania statystyk porownawczych
pobierz_statystyki <- function(wektor_cen, nazwa_zbioru) {
data.frame(
Zbior = nazwa_zbioru,
Min = min(wektor_cen),
Q1 = quantile(wektor_cen, 0.25),
Mediana = median(wektor_cen),
Srednia = mean(wektor_cen),
Q3 = quantile(wektor_cen, 0.75),
Max = max(wektor_cen),
Odchylenie_Std = sd(wektor_cen)
)
}
statystyki_podzialu <- rbind(
pobierz_statystyki(real.train, "Treningowy (70%)"),
pobierz_statystyki(real.test, "Testowy (30%)")
)
cat("--- Porównanie statystyk rozkładu ceny ---\n")## --- Porównanie statystyk rozkładu ceny ---print(statystyki_podzialu)## Zbior Min Q1 Mediana Srednia Q3 Max
## 25% Treningowy (70%) 75000 320000 450000 518500.0 635000 1998000
## 25%1 Testowy (30%) 78000 320000 450000 519955.3 635000 1999950
## Odchylenie_Std
## 25% 285259.2
## 25%1 289769.1Wizualna kontrola gęstości rozkładu cen potwierdza tożsamość obu wydzielonych podzbiorów.
df_wykres_train <- data.frame(Cena = real.train, Zbior = "Treningowy (70%)")
df_wykres_test <- data.frame(Cena = real.test, Zbior = "Testowy (30%)")
df_porownanie <- rbind(df_wykres_train, df_wykres_test)
ggplot(df_porownanie, aes(x = Cena, fill = Zbior)) +
geom_density(alpha = 0.4) +
scale_x_continuous(labels = scales::label_dollar(prefix = "", suffix = " USD"),
limits = c(0, 1500000)) +
scale_fill_manual(values = c("Treningowy (70%)" = "#26A69A", "Testowy (30%)" = "#9C27B0")) +
labs(
title = "Porównanie rozkładu ceny w zbiorze treningowym i testowym",
subtitle = "Weryfikacja zachowania stratyfikacji danych",
x = "Cena nieruchomości",
y = "Gęstość rozkładu",
fill = "Zbiór danych"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 14, hjust = 0.5),
plot.subtitle = element_text(size = 11, hjust = 0.5, color = "gray30"),
legend.position = "bottom"
)
2.3 Skalowanie oraz transformacje predyktorów
W celu eliminacji silnej asymetrii prawostronnej cech metrażowych i cenowych zastosowano transformację logarytmiczną log1p. Dodatkowo, na potrzeby stabilizacji procesu uczenia sieci neuronowych, predyktory numeryczne poddano standaryzacji (Z-score). W celu uniknięcia błędu wycieku danych (data leakage), parametry pozycjonujące i skalujące (średnia i odchylenie standardowe) zostały wyznaczone wyłącznie na próbie treningowej, a następnie zaaplikowane do zbioru testowego.
vars_to_log <- c("price", "sqft_living", "sqft_lot", "sqft_above",
"sqft_basement", "total_sqft", "yr_renovated")
# Zastosowanie transformacji logarytmicznej
housing.train <- housing.train %>% mutate(across(all_of(vars_to_log), log1p))
housing.test <- housing.test %>% mutate(across(all_of(vars_to_log), log1p))
# Aktualizacja wektorów cen o wartości w skali log
real.train <- housing.train$price
real.test <- housing.test$price
# Definicja predyktorów numerycznych do standaryzacji
numeric_predictors <- c("bedrooms", "bathrooms", "sqft_living", "sqft_lot",
"floors", "sqft_above", "sqft_basement", "yr_built",
"yr_renovated", "lat", "long", "house_age",
"total_sqft", "living_ratio")
# Obliczenie profili skalowania wyłącznie na podstawie zbioru treningowego
scaler_profile <- preProcess(housing.train[, numeric_predictors],
method = c("center", "scale"))
# Utworzenie oddzielnych obiektów dedykowanych sieciom neuronowym
housing.train_scaled <- housing.train
housing.test_scaled <- housing.test
# Aplikacja profilu skalującego do zbioru treningowego i testowego
housing.train_scaled[, numeric_predictors] <- predict(scaler_profile, housing.train[, numeric_predictors])
housing.test_scaled[, numeric_predictors] <- predict(scaler_profile, housing.test[, numeric_predictors])
housing.vars.all <- setdiff(names(housing.train), "price")3. Modelowanie i analiza wyników
Użycie struktury kart pozwala na czytelne, niezależne zweryfikowanie każdego z zaimplementowanych algorytmów uczenia nadzorowanego.
# Definicja funkcji metryk przed pierwszym uzyciem w sekcjach modelowania
regresja.metryki <- function(real_log, predicted_log) {
real <- expm1(real_log)
predicted <- pmax(expm1(predicted_log), 0)
text_mse <- mean((real - predicted)^2)
RMSE <- sqrt(text_mse)
MAE <- mean(abs(real - predicted))
MedAE <- median(abs(real - predicted))
MAPE <- mean(abs(real - predicted) / real)
R2 <- cor(predicted, real)^2
wynik <- data.frame(MSE = text_mse, RMSE, MAE, MedAE, MAPE, R2)
return(wynik)
}4. Model Benchmarkowy (OLS)
Klasyczna regresja liniowa estymowana metodą najmniejszych kwadratów stanowi punkt odniesienia dla pozostałych modeli. Model usunął dwie kolumny (house_age oraz is_luxury1) z powodu idealnej współliniowości, co wynika z faktu, że zmienne te zostały sztucznie odtworzone z innych predyktorów (np. roku budowy). Model OLS charakteryzuje się wysoką stabilnością. Współczynnik determinacji wynosi ok. 72% na zbiorze treningowym i 71.5% na zbiorze testowym. Średni błąd procentowy (MAPE) ustabilizował się na poziomie 19.6%. Model liniowy wykazuje jednak ograniczenia strukturalne i nie potrafi wychwycić głębszych interakcji nieliniowych zachodzących na rynku.
# Budujemy klasyczny model regresji liniowej
model.ols <- lm(price ~ . - zipcode, data = housing.train)
# Generowanie predykcji
pred.ols.train <- predict(model.ols, newdata = housing.train)
pred.ols.test <- predict(model.ols, newdata = housing.test)
# Obliczenie metryk
metryki.ols.train <- regresja.metryki(real_log = real.train, predicted_log = pred.ols.train)
metryki.ols.test <- regresja.metryki(real_log = real.test, predicted_log = pred.ols.test)
cat("--- Wyniki regresji OLS na zbiorze testowym ---\n")## --- Wyniki regresji OLS na zbiorze testowym ---print(metryki.ols.test)## MSE RMSE MAE MedAE MAPE R2
## 1 24231482046 155664.6 101959.2 68398.52 0.1964613 0.71568715. Drzewa Decyzyjne
W tej sekcji porównano domyślne drzewo regresyjne z wariantem strojonym metodą walidacji krzyżowej (5-fold CV) pod kątem parametru złożoności. Domyślne drzewo okazało się najsłabszym modelem w zestawieniu (R2=0.628, MAPE = 22.2%), co wynika ze zbyt restrykcyjnych kryteriów stopu algorytmu rpart. Strojenie hiperparametru doprowadziło do wyboru optymalnego CP=0.001. Pozwoliło to na głębszy podział przestrzeni cech i poskutkowało drastyczną poprawą jakości, strojone drzewo osiągnęło R2=0.744 na teście, pokonując tym samym benchmark OLS.
# Drzewo domyślne
tree.default <- rpart(price ~ . - zipcode, data = housing.train, method = "anova")
pred.default.train <- predict(tree.default, newdata = housing.train)
pred.default.test <- predict(tree.default, newdata = housing.test)
metryki.default.train <- regresja.metryki(real_log = real.train, predicted_log = pred.default.train)
metryki.default.test <- regresja.metryki(real_log = real.test, predicted_log = pred.default.test)
# Wizualizacja rpart.plot
layout(matrix(1:1, ncol = 1))
rpart.plot(tree.default, type = 2, extra = 0, fallen.leaves = TRUE, nn = TRUE,
box.palette = "Blues", shadow.col = 0, cex = 0.75,
main = "Domyślne drzewo regresyjne")
# Wizualizacja za pomoca partykit
tree.party <- as.party(tree.default)
plot(tree.party, main = "Drzewo regresyjne z rozkładami cen w liściach",
gp = gpar(fontsize = 7), ip_panel = id_node_inner,
tp_panel = node_boxplot(tree.party, col = "#26A69A", fill = "#26A69A", alpha = 0.5))
# Strojenie drzewa
set.seed(123)
ctrl.tree <- trainControl(method = "cv", number = 5)
grid.tree <- expand.grid(cp = seq(from = 0.001, to = 0.05, by = 0.002))
tree.tuned <- train(price ~ . - zipcode, data = housing.train, method = "rpart",
trControl = ctrl.tree, tuneGrid = grid.tree, metric = "RMSE")
pred.tree.tuned.train <- predict(tree.tuned, newdata = housing.train)
pred.tree.tuned.test <- predict(tree.tuned, newdata = housing.test)
metryki.tree.tuned.train <- regresja.metryki(real_log = real.train, predicted_log = pred.tree.tuned.train)
metryki.tree.tuned.test <- regresja.metryki(real_log = real.test, predicted_log = pred.tree.tuned.test)
cat("--- Wyniki strojonego drzewa na zbiorze testowym ---\n")## --- Wyniki strojonego drzewa na zbiorze testowym ---print(metryki.tree.tuned.test)## MSE RMSE MAE MedAE MAPE R2
## 1 21853827364 147830.4 92444.39 56535.6 0.177291 0.7442266. Las Losowy (Random Forest)
Metody zespołowe typu Bagging i Random Forest redukują wariancję pojedynczych drzew poprzez uśrednianie predykcji. Domyślny las losowy (100 drzew) osiągnął bardzo wysokie dopasowanie do danych testowych (R2=0.856, MAPE = 13.18%). Model strojony za pomocą szybkiego silnika ranger wskazał jako optymalne parametry: mtry = 9 oraz min.node.size = 5. Zgodnie z oczekiwaniami teoretycznymi, algorytm silnie dopasowuje się do próby treningowej (R2=0.969), wykazując cechy przeuczenia. Mimo to, jego zdolność generalizacji na zbiorze testowym pozostaje istotna (R2=0.852). Wykres ważności cech wskazuje, że kluczowy wpływ na cenę mają: wielkość powierzchni mieszkalnej (sqft_living), standard wykończenia (grade) oraz lokalizacja geograficzna (lat).
# Las domyślny
set.seed(123)
rf.default <- randomForest(price ~ . - zipcode, data = housing.train, ntree = 100, importance = TRUE)
pred.rf.default.train <- predict(rf.default, newdata = housing.train)
pred.rf.default.test <- predict(rf.default, newdata = housing.test)
metryki.rf.default.train <- regresja.metryki(real_log = real.train, predicted_log = pred.rf.default.train)
metryki.rf.default.test <- regresja.metryki(real_log = real.test, predicted_log = pred.rf.default.test)
# Las strojony
set.seed(123)
ctrl.rf <- trainControl(method = "cv", number = 5)
grid.rf <- expand.grid(mtry = c(3, 6, 9), min.node.size = c(5, 15), splitrule = "variance")
rf.tuned <- train(price ~ . - zipcode, data = housing.train, method = "ranger",
trControl = ctrl.rf, tuneGrid = grid.rf, num.trees = 100, importance = "permutation")
pred.rf.tuned.train <- predict(rf.tuned, newdata = housing.train)
pred.rf.tuned.test <- predict(rf.tuned, newdata = housing.test)
metryki.rf.tuned.train <- regresja.metryki(real_log = real.train, predicted_log = pred.rf.tuned.train)
metryki.rf.tuned.test <- regresja.metryki(real_log = real.test, predicted_log = pred.rf.tuned.test)
# Wykres ważności zmiennych
importance.rf <- varImp(rf.tuned)$importance
importance.df <- data.frame(Zmienna = rownames(importance.rf), Waznosc = importance.rf$Overall) %>% arrange(desc(Waznosc))
ggplot(importance.df, aes(x = reorder(Zmienna, Waznosc), y = Waznosc)) +
geom_bar(stat = "identity", fill = "#26A69A", alpha = 0.8) +
coord_flip() +
labs(title = "Ważność zmiennych w optymalnym modelu Random Forest", x = "Zmienne predykcyjne", y = "Wartość ważności") +
theme_minimal()
# Wykres diagnostyczny Real vs Predicted
df_pred_vs_real <- data.frame(Rzeczywiste = expm1(real.test), Przewidywane = expm1(pred.rf.tuned.test))
ggplot(df_pred_vs_real, aes(x = Rzeczywiste, y = Przewidywane)) +
geom_point(alpha = 0.4, color = "#26A69A") +
geom_abline(intercept = 0, slope = 1, color = "#FF7043", linetype = "dashed", linewidth = 1) +
scale_x_continuous(labels = scales::label_dollar(prefix = "", suffix = " USD"), limits = c(0, 1500000)) +
scale_y_continuous(labels = scales::label_dollar(prefix = "", suffix = " USD"), limits = c(0, 1500000)) +
labs(title = "Wartości rzeczywiste vs przewidywane (Random Forest)", x = "Realna cena", y = "Cena przewidziana") +
theme_minimal()
7. Sieci Neuronowe
Sztuczne sieci neuronowe zaimplementowane za pomocą pakietu neuralnet wymagały pełnego przekształcenia wejściowych faktorów na postać zero-jedynkową (One-Hot Encoding). Przetestowano dwa warianty architektury: sieć jednowarstwową z 6 neuronami ukrytymi c(6) oraz sieć głębszą, dwuwarstwową c(4, 2). Wariant jednowarstwowy NN_c6 okazał się najlepszym modelem pod względem błędu globalnego. Potrzebował około 30 minut na zbieżność algorytmu na pełnym zbiorze, ale wygenerował najniższy błąd średniokwadratowy na próbie ukrytej (RMSE = 112 577 USD). Co niezwykle istotne z punktu widzenia metodologicznego, model ten wykazuje minimalny, niemal zerowy rozjazd między metrykami treningowymi (R2=0.865) a testowymi (R 2=0.850). Świadczy to o bardzo dobrej generalizacji wiedzy i optymalnym doborze struktur sieci, bez zjawiska szkodliwego przeuczenia. Wariant głębszy NN_c4_2 okazał się zbyt wąski w kolejnych warstwach i odnotował słabsze wyniki (R 2=0.833).
# Przygotowanie danych (One-Hot Encoding)
X_train_matrix <- model.matrix(price ~ . - zipcode, data = housing.train_scaled)
X_test_matrix <- model.matrix(price ~ . - zipcode, data = housing.test_scaled)
X_train_nn <- as.data.frame(X_train_matrix)
X_test_nn <- as.data.frame(X_test_matrix)
names(X_train_nn) <- make.names(names(X_train_nn))
names(X_test_nn) <- make.names(names(X_test_nn))
X_train_nn$price <- housing.train_scaled$price
X_test_nn$price <- housing.test_scaled$price
nn_predictors <- setdiff(names(X_train_nn), c("X.Intercept.", "price"))
# Ładowanie modeli z sesji w celu szybkiego renderowania raportu
nn_housing <- readRDS("wytrenowane_modele/housing.nn.c6.rds")
nn_housing_2 <- readRDS("wytrenowane_modele/housing.nn.c4.2.rds")
# Wyciągnięcie wyników na zbiorze testowym
pred_nn_test_log <- compute(nn_housing, X_test_nn[, nn_predictors])$net.result
metryki.nn.test <- regresja.metryki(real_log = real.test, predicted_log = as.numeric(pred_nn_test_log))
pred_nn2_test_log <- compute(nn_housing_2, X_test_nn[, nn_predictors])$net.result
metryki.nn2.test <- regresja.metryki(real_log = real.test, predicted_log = as.numeric(pred_nn2_test_log))
# Wizualizacja sieci c(6)
plotnet(nn_housing, circle_col = list("lightyellow", "lightblue"), bord_col = "grey40",
pos_col = "darkgreen", neg_col = "red3", circle_cex = 2, cex_val = 0.5, max_sp = TRUE,
main = "Topologia sieci neuronowej c(6)")
8. Podsumowanie i porównanie modeli
8.1 Zbiorcze zestawienie metryk
Poniższe tabele prezentują kompletne zestawienie wszystkich zbudowanych modeli, odtworzone bezpośrednio z uzyskanych wyników empirycznych.
8.1.1 Wyniki na zbiorze Treningowym
# Reczne odtworzenie tabeli zbiorczej
tabela_train <- data.frame(
Model = c("OLS", "Tree_Default", "Tree_Tuned", "RF_Default", "RF_Tuned", "NN_c6", "NN_c4_2"),
MSE = c(23017421073, 29644177967, 18984669054, 2752151890, 3008391760, 11102941971, 13226030478),
RMSE = c(151714.93, 172174.85, 137784.87, 52460.96, 54848.81, 105370.50, 115004.48),
MAE = c(100299.77, 111602.74, 87894.60, 30094.76, 31453.75, 67891.60, 72612.29),
MedAE = c(67363.04, 71790.04, 55641.31, 17221.77, 17829.84, 42864.26, 45199.86),
MAPE = c(0.1957, 0.2161, 0.1707, 0.0573, 0.0597, 0.1335, 0.1410),
R2 = c(0.7200, 0.6420, 0.7695, 0.9717, 0.9690, 0.8647, 0.8387)
)
knitr::kable(tabela_train, digits = 4, format = "markdown")| Model | MSE | RMSE | MAE | MedAE | MAPE | R2 |
|---|---|---|---|---|---|---|
| OLS | 23017421073 | 151714.93 | 100299.77 | 67363.04 | 0.1957 | 0.7200 |
| Tree_Default | 29644177967 | 172174.85 | 111602.74 | 71790.04 | 0.2161 | 0.6420 |
| Tree_Tuned | 18984669054 | 137784.87 | 87894.60 | 55641.31 | 0.1707 | 0.7695 |
| RF_Default | 2752151890 | 52460.96 | 30094.76 | 17221.77 | 0.0573 | 0.9717 |
| RF_Tuned | 3008391760 | 54848.81 | 31453.75 | 17829.84 | 0.0597 | 0.9690 |
| NN_c6 | 11102941971 | 105370.50 | 67891.60 | 42864.26 | 0.1335 | 0.8647 |
| NN_c4_2 | 13226030478 | 115004.48 | 72612.29 | 45199.86 | 0.1410 | 0.8387 |
8.1.2 Wyniki na zbiorze Testowym
tabela_test <- data.frame(
Model = c("NN_c6", "RF_Default", "RF_Tuned", "NN_c4_2", "Tree_Tuned", "OLS", "Tree_Default"),
MSE = c(12673581852, 12905411967, 13307955582, 14110750344, 21853827364, 24231482046, 31920335784),
RMSE = c(112577.0, 113602.0, 115360.1, 118788.7, 147830.4, 155664.6, 178662.6),
MAE = c(70609.89, 68385.35, 69424.56, 74371.99, 92444.39, 101959.17, 115402.87),
MedAE = c(44658.62, 40261.91, 41046.83, 46019.75, 56535.60, 68398.52, 73211.88),
MAPE = c(0.1370, 0.1318, 0.1336, 0.1434, 0.1773, 0.1965, 0.2219),
R2 = c(0.8503, 0.8561, 0.8519, 0.8336, 0.7442, 0.7157, 0.6280)
)
knitr::kable(tabela_test, digits = 4, format = "markdown")| Model | MSE | RMSE | MAE | MedAE | MAPE | R2 |
|---|---|---|---|---|---|---|
| NN_c6 | 12673581852 | 112577.0 | 70609.89 | 44658.62 | 0.1370 | 0.8503 |
| RF_Default | 12905411967 | 113602.0 | 68385.35 | 40261.91 | 0.1318 | 0.8561 |
| RF_Tuned | 13307955582 | 115360.1 | 69424.56 | 41046.83 | 0.1336 | 0.8519 |
| NN_c4_2 | 14110750344 | 118788.7 | 74371.99 | 46019.75 | 0.1434 | 0.8336 |
| Tree_Tuned | 21853827364 | 147830.4 | 92444.39 | 56535.60 | 0.1773 | 0.7442 |
| OLS | 24231482046 | 155664.6 | 101959.17 | 68398.52 | 0.1965 | 0.7157 |
| Tree_Default | 31920335784 | 178662.6 | 115402.87 | 73211.88 | 0.2219 | 0.6280 |
8.2 Wnioski końcowe i rekomendacje biznesowe
Model sieci neuronowej o architekturze NN_c6 okazał się najlepszym narzędziem estymacyjnym na zbiorze ukrytym, osiągając najniższy błąd globalny (RMSE = 112 577 USD) oraz wysoki poziom wyjaśnionej wariancji (R2=85.03%). Drugie miejsce zajął domyślny las losowy (RF_Default) z minimalnie lepszym współczynnikiem determinacji (R2=85.61%).
Wszystkie modele zaawansowane (z wyjątkiem sztywnego, domyślnego drzewa rpart) bezproblemowo pobiły klasyczny benchmark ekonometryczny OLS (R2=0.7157). Pokazuje to, że struktura zależności cenowych na rynku nieruchomości opiera się na silnych, nieliniowych interakcjach cech przestrzennych i standardu budynków.
Modele lasów losowych osiągnęły wysokie dopasowanie na zbiorze treningowym (R2 powyżej 96.8%), jednak na zbiorze testowym zanotowały nieznaczny spadek dopasowania, ustępując miejsca stabilnej i zrównoważonej sieci neuronowej, która wykazuje minimalne zróżnicowanie między treningiem a testem.
Osiągnięcie średniego procentowego błędu prognozy na poziomie 13.1% - 13.7% przez modele zespołowe i sieciowe jest wynikiem w pełni satysfakcjonującym w standardach rynkowych. Opracowany system może z powodzeniem wspierać procesy masowej, automatycznej wyceny nieruchomości (AVM) dla bankowości oraz agencji pośrednictwa.
# Koncowe zabezpieczenie danych na dysku
save(housing.train, housing.test, housing.train_scaled, housing.test_scaled,
housing.vars.all, real.train, real.test, file = "housing_prepared.RData")