Przewidywanie głosów za pomocą modeli liniowych

Ta część pokazuje, jak pracować z modelami statystycznymi przy użyciu R. Pokazuje, jak sprawdzić założenia dotyczące danych, określić modele liniowe, tworzyć prognozy i mierzyć dokładność predykcyjną. Pokazuje również, jak programowo znaleźć dobre modele, aby uniknąć ręcznego wykonywania analiz, co może potencjalnie zaoszczędzić dużo czasu. Pod koniec będziemy pracować z różnymi narzędziami ilościowymi, które są obecnie używane w wielu obszarach biznesowych i badawczych. Pakiety użyte tu  są tymi samymi, co w poprzedniej części. Podobnie jak  poprzednio, skupimy się tutaj na programowej automatyzacji analizy, a nie na dogłębnym zrozumieniu technik statystycznych zastosowanych w tej części. Omówimy następujące kwestie:

* Dzielenie danych na zbiory szkoleniowe i testowe

* Tworzenie modeli regresji liniowej służących do prognozowania

* Sprawdzanie założeń modelu różnymi technikami

* Pomiar dokładności predykcyjnej dla danych liczbowych i jakościowych

* Programowe wyszukiwanie najlepszego możliwego modelu

Wymagane pakiety

ggplot2 –  Wysokiej jakości wykresy

corrplot –  Wykresy korelacji

progres – Pokaż postęp iteracji

Programowanie poprzez wizualizację pełnego obrazu

Teraz będziemy pracować z podejściem odgórnym, co oznacza, że ​​najpierw zaczniemy od kodu abstrakcyjnego i stopniowo przejdziemy do szczegółów implementacji. Generalnie uważam, że to podejście jest bardziej wydajne, gdy masz jasne wyobrażenie o tym, co chcesz zrobić. W naszym przypadku zaczniemy pracując z plikiem main.R. Pierwszą rzeczą, na którą należy zwrócić uwagę, jest to, że użyjemy funkcji proc.time() dwukrotnie, raz na początku i raz na końcu, i użyjemy różnicy między tymi dwiema wartościami, aby zmierzyć, ile czasu wykonanie całego kodu zajęło. Drugą rzeczą, na którą należy zwrócić uwagę, jest to, że funkcja empty_directories() zapewnia istnienie każdego z określonych katalogów i usuwa wszystkie zawarte w nich pliki. Używamy go do czyszczenia naszych katalogów na początku każdego wykonania, aby upewnić się, że mamy najnowsze pliki i tylko pliki utworzone w ostatnim uruchomieniu. Rzeczywisty kod jest pokazany poniżej i po prostu iteruje przez każdy z przekazanych katalogów, rekurencyjnie usuwa wszystkie pliki wewnątrz za pomocą funkcji unlink() i upewnia się, że katalog istnieje z funkcji dir.create(). Unika wyświetlania jakichkolwiek ostrzeżeń związanych z katalogiem już istniejącym, co nie stanowi problemu w naszym przypadku, poprzez użycie parametru showWarnings = FALSE.

empty_directories() <- function(directories) {

for (directory in directories) {

unlink(directory, recursive = TRUE)

dir.create(driectory, showWarnings = FALSE)

}

}

Od początku R, używamy funkcji print_section() i  empty_directories() do wyświetlania nagłówków i usuwania zawartości katalogu (w celu odtworzenia wyników za każdym razem, gdy uruchamiamy funkcję z pustymi katalogami), a użyjemy mechanizm pokazany za pomocą proc.time() do pomiaru czasu wykonania. Teraz, gdy poprzednie dwa punkty są już na uboczu, przejdziemy do pokazania pełnej zawartości pliku main.R.

start_time <- proc_time()

source(„./function.R”)

empty_directories(c (

„./results/oroginal/”,

„./results/adjusted/’,

„./results/original/scatter_plots/”

))

data <- prepare_data(„./data_brexit_referendum.csv, complete_cases = TRUE)

data_adjusted  <- adjust_data(data)

numerical_variables <- get_numerical_variable_names(data)

numerical_variables_adj <- get_numerical_variable_names(data_adjusted)

print(„Working on summaries …”)

full_summary(data, save_to = „./reults/original/summary_text.txt”)

numerical_summary)

data,

numerical_variables = numerical_variables,

save_to = „./results/original/summary_numerical.csv”

}

print(„Working on histograms …”)

plot_percentage)

data,

variable = „RegionName”

save_to = „./results/origibnal/vote_percentage_by-egion.png”

}

print („Working on matrix scatter plots…”)

matrix_scatter_plots (

data_adjusted,

nuerical_variables = numerical_variables_adj,

save_to = „./results/adjusted/matrix_scatter_plots.png”

)

print(„Working on scatter plots …”)

plot_scatter_plot (

data;

var_x = „RegionName”,

var_y = „Proportion”,

var_color = „White”,

regression = TRUE,

save_to = „./results/origina;/reginname_vs_proportion_vs_white.png”

)

all_scatter_plots(

data,

numerical_variables = numerical_variables,

save_to = „./results/original/scatter_plots/”

)

print („Working on correlations…”)

correlations_plot(

data,

numerical_variables = numerical_variables,

save_to = „./results/original/correlations.png”

)

print(„Working on principal components…”)

principal_components(

data_adjusted,

numerical_variables – numerical_variables_adj,

save_to = „./results/adjusted/principal_components”

)

end_time <- proc.time()

time_taken <- end_time – start.time

print(paste(„Time taken: „ , taken[1]))

print („Done.”)

Jak widać, mając tylko ten plik, masz duży obraz analizy i jesteś w stanie odtworzyć swoją analizę, uruchamiając pojedynczy plik, zapisując wyniki na dysku (zwróć uwagę na argumenty save_to,), i zmierzyć czas potrzebny do wykonania pełnej analizy. Z naszej listy celów ogólnych, ten kod spełnia cele od pierwszego do czwartego. Osiągnięcie celu piątego i szóstego zostanie osiągnięte poprzez pracę nad plikiem function.R który zawiera wiele małych funkcji. Posiadanie tego pliku main.R daje nam mapę tego, co należy zaprogramować, i chociaż w tej chwili to nie zadziała, ponieważ funkcje, których używa, jeszcze nie istnieją, zanim zakończymy ich programowanie, plik ten będzie nie wymagają żadnych zmian i przyniosą pożądane rezultaty. Ze względu na ograniczenia miejsca nie przyjrzymy się implementacji wszystkich funkcji w pliku main.R, tylko te reprezentatywne: prepare.data(), plot_scatter_plot() i all_scatter_plots(). Inne funkcje wykorzystują podobne techniki do hermetyzacji odpowiedniego kodu. Zaczynamy od prepare_data(). Ta funkcja jest abstrakcyjna i wykorzystuje cztery różne konkretne funkcje do wykonania swojej pracy: read_csv(), clean_data(), transform_data() i, jeśli jest to wymagane, complete.cases(). Pierwsza funkcja, mianowicie read.csv(), odbiera ścieżkę do pliku CSV, z którego dane są odczytywane i ładuje do obiektu ramki danych o nazwie data w tym przypadku. Czwarta funkcja, widzieliście wcześniej , więc nie będziemy jej tutaj wyjaśniać. Funkcje dwa i trzy zostały stworzone przez nas i wyjaśnimy je. Zauważ, że main.R nie wie o tym, jak przygotowywane są dane, prosi tylko o przygotowanie danych i deleguje zadanie do funkcji abstrakcyjnej prepare_data()

prepare_data <- function(path, complete_cases = TRUE) {

data <- read.csv(path)

data <- clean_data(data)

data <- transform_data(data)

if (complete_cases) {

data <- data[complete.cases(data), ]

}

return(data)

}

Funkcja clean_data() po prostu hermetyzuje na razie ponowne kodowanie -1 dla NA. Jeśli nasza procedura czyszczenia nagle stała się bardziej złożona (na przykład nowe źródła danych wymagające częstszego czyszczenia lub uświadomienia sobie, że coś przeoczyliśmy i musimy dodać to do procedury czyszczenia), dodalibyśmy te zmiany do tej funkcji i nie musielibyśmy modyfikować niczego innego w pozostałej części naszego kodu. Oto niektóre z zalet hermetyzacji kodu w funkcje, które przekazują zamiar i izolują to, co należy zrobić w małych krokach:

clean_data <- function(data) {

data[data$Leave == -1, „Leave”] <- NA

return(data)

}

Aby przekształcić nasze dane, dodając dodatkowe zmienne Proportion i Vote , i ponownie oznacz plik nazwy regionów, używamy następującej funkcji:

transform_data <- funtion(data) {

data$Proportion <- data$Leave / data$NVotes

data$Vote <- ifelse(data$Proportion > 0,5 , „Leave”, „Remain”)

data$RegionName <- as.character(data$RegionName)

data[data$RegionName == „London”, „RegionName”] <- „L”

data[data$RegionName == „North West”, „RegionName”] <- „NW”

data[data$RegionName == „South West”, „RegionName”] <- „SW”

data[data$RegionName == „South East”, „RegionName”] <- „SE”

data[data$RegionName == „East Middlands”, „RegionName”] <- „EM”

data[data$RegionName == „West Midlands”, „RegionName”] <- „WM”

data[data$RegionName == „East of England” , „RegionName”] <- „EE”

data[data$RegionName == „Yorkshire and The Humer”, „RegionName”] <- „Y”

return(data)

}

Wszystkie te linie kodu, które widzieliście wcześniej. Wszystko, co robimy, to hermetyzowanie ich w funkcje, które komunikują intencje i pozwalają nam znaleźć miejsce, w którym zachodzą określone procedury, abyśmy mogli je znaleźć i łatwo zmienić, jeśli zajdzie taka potrzeba później.

Teraz przyjrzymy się plot_scatter_plot(). Ta funkcja znajduje się pomiędzy funkcją abstrakcyjną a konkretną. Użyjemy go bezpośrednio w naszym pliku main.R, ale użyjemy go również w innych funkcjach w pliku function.R. Wiemy, że przez większość czasu będziemy używać Proportion jako zmiennej koloru, więc dodajemy ją jako wartość domyślną, ale pozwalamy użytkownikowi całkowicie usunąć kolor, sprawdzając, czy argument został wysłany jako FALSE, a ponieważ będziemy używać tej samej funkcji do tworzenia wykresów, które przypominają wszystkie wykresy rozrzutu, utworzonych do tego momentu, uczynimy linię regresji opcjonalną. Zwróć uwagę, że w przypadku pierwszych wykresów oś x jest zmienną ciągłą, ale w przypadku drugiego wykresu jest to zmienna jakościowa (czynnikowa). Ten rodzaj elastyczności jest bardzo potężny i jest dostępny dla nas ze względu na zdolność ggplot2 do dostosowania się do tych zmian. Formalnie nazywa się to polimorfizmem. Wreszcie, zamiast zakładać, że użytkownik zawsze będzie chciał zapisać wynikowy wykres na dysku, ustawiamy argument save_to jako opcjonalny, podając dla niego pusty ciąg. W razie potrzeby sprawdzamy, czy ten ciąg jest pusty za pomocą not_empty(), a jeśli nie jest pusty, konfigurujemy mechanizm zapisywania PNG.

plot_scatter_plot <- function(datal

var_x,

var_y,

var_color = „Proportion”,

regression = FALSE;

save_to = „ „ ) {

if (var_color) {

plot <- ggplot(data, aes_string(x = var_x, y = var_y, color = var_color))

} else {

plot <- ggplot(data, aes_string(x = var_x, y = var_y))

}

plot <- plot + scale_color_viridis()

plot <- plot + geom_point()

if (regression) {

plot <- plot + stat_smooth(method = „lm”, col = „grey”, se = FALSE)

}

if (not_empty(save_to)) png(save_to)

print(plot)

if (not_empty (save_to)) dev.off()

}

Teraz przyjrzymy się all_scatter_plots(). Ta funkcja jest funkcją abstrakcyjną, która ukrywa przed wiedzą użytkownika nazwę funkcji, która iteracyjnie tworzy wykresy, dogodnie nazwana create_graphs_iteratively() i funkcja graficzna, funkcja plot_scatter_plot(), którą widzieliśmy wcześniej. W przypadku, gdy chcemy ulepszyć mechanizm iteracyjny lub funkcję graficzną, możemy to zrobić bez konieczności wprowadzania zmian od osób, które używają naszego kodu, ponieważ ta wiedza jest tutaj zamknięta. Funkcja create_graphs_iteratively() jest taka sama, jak wcześniej, z wyjątkiem kodu paskowego postępu. Pakiet progress udostępnia funkcję progres_bar$new(), która tworzy pasek postępu w terminalu podczas wykonywania iteracyjnego procesu, dzięki czemu widzimy jaki procent procesu został ukończony i dowiedz się, ile czasu pozostało. Zwróć uwagę na zmianę argumentu save_to z funkcji plot_scatter_plot() i all_scatter_plot(). W pierwszym przypadku jest to nazwa pliku; w drugim – nazwa katalogu. Różnica jest niewielka, ale ważna. Nieostrożny czytelnik może tego nie zauważyć i może być przyczyną nieporozumień. Funkcja plot_scatter_plot() tworzy pojedynczy wykres i otrzymuje w ten sposób nazwę pliku. Jednak all_scatter_plots() wyprodukuje, używając plot_scatter_plot(),  dużo wykresów, więc musi wiedzieć, gdzie trzeba je wszystkie zapisać, dynamicznie tworzyć ostateczne nazwy obrazów i wysyłać je pojedynczo do plot_scatter_plot().

Wreszcie, ponieważ chcemy, aby regresja została uwzględniona na tych wykresach, po prostu wysyłamy parametr regression = TREU:

all_scatter-plots <- function(data, numerical_variables, save_to = „ „ ) {

create_graphs_iteratively(data, numerical_variables, plot_scatter_plot, save_to)

}

create_graphs_iteratively <- function(data,

numerical_variables,

plot_function,

save_to = „ „) {

numerical_variables[[„Proportion”]] <- FALSE

variables <- names(numerical_Variables[numerical_variables == TRUE])

n_variables <- (length(variables) – 1)

progres_bar <- progres_bar$new(

format = „Progress [:bar] :percent ETA :  :ets”,

total = n_Variables

)

for (i in 1:n_variables) {

progres_bar$tick()

for (j in (i+1) : length(variables)) {

image_name <- paste(

save_to,

variables[i], „_”,

variables[j], „.png”,

sep= „”

)

plot_function(

data,

var_x = variables[i].

var_y = variables[j],

save_to = „image_name,

regression = TRUE

}

}

}

}

Zrozumienie podstaw kodu wysokiej jakości

Kod, który jest modułowy, elastyczny i którego zależności są dobrze zarządzane, jest uważany za wysoce spójny i luźno powiązany. Terminy te są najczęściej używane w środowiskach zorientowanych obiektowo, ale ogólnie mają zastosowanie do każdego systemu. Wysoce spójne oznacza, że ​​rzeczy, które mają być razem, są. Luźno powiązane oznacza, że ​​rzeczy, które nie powinny być razem, nie są. Poniższy rysunek przedstawia te cechy, w których każdy z okręgów może być funkcją lub ogólnie obiektem. To są podstawy zarządzania zależnościami. Wiele książek poświęconych tym zagadnieniom zostało i nadal jest publikowanych.

Najważniejsze zasady dotyczące kodu wysokiej jakości to:

  1. Zrób rzeczy małe i skupione na jednej odpowiedzialności.
  2. Niech konkrety zależą od abstrakcji (nie odwrotnie).
  3. Twórz rzeczy, które są wysoce spójne i luźno powiązane.

Zaczynamy od utworzenia dwóch plików: functions.R i main.R. Plik function.R zawiera funkcje wysokiego poziomu (głównie wywoływane z pliku main.R), a także funkcje niskiego poziomu (używane w innych funkcjach). Czytając plik main.R, powinniśmy mieć jasne wyobrażenie o tym, co robi analiza (taki jest cel funkcji wysokiego poziomu), a wykonanie jej powinno odtworzyć naszą analizę dla wszelkich danych, które pasują do naszych podstawowe założenia (w tym przykładzie są to głównie struktury danych). Powinniśmy zawsze utrzymywać powiązany kod na tym samym poziomie abstrakcji. Oznacza to, że nie chcemy programować rzeczy na ogólnym poziomie i implementować ich z mieszanymi szczegółami, a rozdzielanie naszego kodu na main.R i function.R jest pierwszy krok w tym kierunku. Ponadto żaden kod w pliku main.R nie powinien zależeć od szczegółów implementacji. To czyni go modułowym w tym sensie, że jeśli chcemy zmienić sposób implementacji czegoś, możemy to zrobić bez konieczności zmiany kodu wysokiego poziomu. Jednak sposób, w jaki implementujemy rzeczy, zależy od tego, co chcemy, aby analiza ostatecznie zrobiła, co oznacza, że ​​konkretne implementacje powinny zależeć od implementacji abstrakcyjnych, które z kolei zależą od celu naszej analizy (określonego jako kod w pliku main.R). Kiedy przenosimy wiedzę z jednego zestawu kodu do innego, generujemy zależność, ponieważ kod, który wie o innym kodzie, zależy od tego, czy działa poprawnie. Chcemy jak najbardziej unikać tych zależności, a co najważniejsze – kierować ich kierunkiem. Jak stwierdzono wcześniej, streszczenie nie powinno zależeć od konkretu, ani inaczej mówiąc, konkret powinien zależeć od abstraktu. Ponieważ analiza (main.R) jest po stronie abstrakcyjnej nie powinno zależeć od szczegółów implementacji konkretnych funkcji. Ale jak można przeprowadzić naszą analizę bez znajomości funkcji, które ją wdrażają? Cóż, nie może. Dlatego potrzebujemy pośrednika, funkcji abstrakcyjnych. Funkcje te mają na celu main.R dostarczenie stabilnej wiedzy i zagwarantowanie, że analiza, której poszukuje, zostanie wykonana, a także eliminują zależność implementacji od szczegółów implementacji poprzez samodzielne zarządzanie tą wiedzą. Może się to wydawać zawiłym sposobem pracy i trudną koncepcją do zrozumienia, ale kiedy to zrobisz, przekonasz się, że jest to bardzo proste i będziesz w stanie stworzyć kod, który można podłączyć, co jest dużym wzrostem wydajności. Możesz zajrzeć do wcześniej wspomnianych książek, aby uzyskać głębszy sens tych koncepcji.

Planowanie przed programowaniem

Często ludzie rozpoczynają programowanie, zanim mają ogólne pojęcie o tym, co chcą osiągnąć. Jeśli jesteś doświadczonym programistą, może to być dobry sposób na zrozumienie problemu, ponieważ masz już rozwiniętą intuicję i prawdopodobnie i tak wyrzucisz kilka pierwszych prób. Jeśli jednak jesteś początkującym programistą, radzę jasno i jasno sformułować swoje cele przed napisaniem jakiegokolwiek kodu (wprowadzenie ich na piśmie może pomóc). Pomoże Ci to w podejmowaniu lepszych decyzji, zadając sobie pytanie, w jaki sposób określony sposób działania wpłynie na Twoje cele. Dlatego zanim cokolwiek ustalimy, musimy zrozumieć i jasno określić nasze ogólne cele:

  1. Szybko zrozum ogólny obraz analizy.
  2. Automatycznie odtworzyć naszą analizę, wykonując pojedynczy plik.
  3. Zapisz wszystkie wynikowe obiekty, tekst i obrazy do analizy.
  4. Zmierz ilość czasu potrzebnego do wykonania pełnej analizy.
  5. Pracując nad procesami iteracyjnymi, poznaj procent ukończenia.
  6. Umieć łatwo znaleźć i zmienić każdą część analizy.

Aby spełnić te ogólne cele, musimy opracować kod modułowy z dobrze zarządzanymi zależnościami, które są elastyczne (łatwe do zmiany) i przyjazne dla efektów ubocznych (zapisywanie obiektów, tekstów i obrazów). Nawet jeśli twoje wyraźne cele tego nie wymagają, powinieneś nabrać zwyczaju programowania w ten sposób, nawet jeśli robisz tylko analizę danych.

Wszystko razem w wysokiej jakości kod

Teraz, gdy mamy już podstawy dotyczące analizowania danych za pomocą statystyk opisowych, zamierzamy ulepszyć strukturę i elastyczność naszego kodu, dzieląc go na funkcje. Chociaż jest to powszechna wiedza wśród wydajnych programistów, nie jest to powszechna praktyka wśród analityków danych. Wielu analityków danych po prostu wkleiłoby kod, który opracowaliśmy razem, w takiej postaci, w jakiej jest, do jednego pliku i uruchamiałoby go za każdym razem, gdy chcieli przeprowadzić analizę. Nie będziemy dodawać nowych funkcji do analizy. Wszystko, co zrobimy, to zmienić kolejność kodu w funkcje, aby hermetyzować ich wewnętrzne działanie i przekazywać intencje za pomocą nazw funkcji (to znacznie zmniejsza potrzebę komentarzy). Skoncentrujemy się na tworzeniu kodu wysokiej jakości, który jest łatwy do odczytania, ponownego użycia, modyfikacji i naprawy (w przypadku błędów). Sposób, w jaki to robimy, jest kwestią stylu, a różne sposoby aranżacji kodu pasują do różnych kontekstów. Metoda, z którą będziemy tutaj pracować, dobrze mi służyła w różnych sytuacjach, ale może nie być najlepsza dla Ciebie. Jeśli nie odpowiada Twoim potrzebom, możesz to zmienić. Niezależnie od preferowanego stylu, inwestycja w tworzenie nawyku ciągłego tworzenia kodu wysokiej jakości sprawi, że na dłuższą metę staniesz się bardziej wydajnym programistą i nadejdzie moment, w którym nie będziesz chciał programować nieefektywnie.

Budowanie nowych zmiennych za pomocą głównych komponentów

Analiza głównych składowych (PCA) to technika redukcji wymiarowości, która jest szeroko stosowana w analizie danych, gdy istnieje wiele zmiennych numerycznych, z których niektóre mogą być skorelowane, i chcielibyśmy zmniejszyć liczbę wymiarów wymaganych do zrozumienia danych. Pomocne może być zrozumienie danych, ponieważ myślenie w więcej niż trzech wymiarach może być problematyczne, a także przyspieszenie algorytmów wymagających dużej mocy obliczeniowej, zwłaszcza w przypadku dużej liczby zmiennych. Dzięki PCA możemy wyodrębnić większość informacji tylko do jednej lub dwóch zmiennych skonstruowanych w bardzo specyficzny sposób, tak aby wychwytywały większość wariancji, a jednocześnie miały dodatkową korzyść w postaci braku korelacji między nimi przez konstrukcję.

Pierwszy główny składnik to liniowa kombinacja oryginalnych zmiennych, która oddaje maksymalną wariancję (informacje) w zbiorze danych. Żaden inny składnik nie może mieć większej zmienności niż pierwszy składnik główny. Następnie druga składowa główna jest ortogonalna do pierwszej i jest obliczana w taki sposób, aby uchwycić maksymalną wariancję pozostałą w danych. I tak dalej. Fakt, że wszystkie zmienne są kombinacjami liniowymi, które są między sobą ortogonalne, jest kluczem do ich nieskorelowania między sobą.

Wystarczy mówić o statystykach; przejdźmy do programowania! Wykonując PCA w R mamy wiele funkcji, które mogą wykonać to zadanie. Aby wspomnieć o niektórych z nich, mamy prcomp() i princom() z pakietu stats (wbudowanego), PCA() z pakietu FactoMineR , dudi.pca() z pakietu ade4 i acp() z pakietu amap. W naszym przypadku użyjemy funkcji prcomp() wbudowanej w R. Aby wykonać PCA, użyjemy skorygowanych danych z poprzedniej sekcji. Najpierw usuwamy zmienne numeryczne, które są skorelowane z Proportion. Następnie wykonujemy PCA, wysyłając dane liczbowe do funkcji prcomp(), a także niektóre parametry normalizacji. center . = TRUE odejmie od siebie średnią każdej zmiennej, oraz scale. = TRUE sprawi, że wariancja każdej zmiennej będzie jednolita, skutecznie normalizując dane. Normalizacja danych jest bardzo ważna podczas wykonywania PCA, ponieważ jest to metoda wrażliwa na skale:

numerical_variables_adjusted[[„Votes”]] <- FALSE

numerical_varaibles_adjusted[[„Leave”]] <- FALSE

data_numeircal_adjusted <- data_adjusted[, numerical_variables_adjusted]

pca <- promp(data_numerical_adjusted, center = TRUE, scale. = TRUE)

pca

#> Odchylenia standardowe (1, .., p = 21):

#> [1] 2,93919 2,42551 1,25860 1,13300 1,00800 0,94112 0,71392 0,57613

#> [9] 0,54047 0,44767 0,37701 0,30166 0,21211 0,17316 0,13759 0,11474

#> [17] 0,10843 0,09797 0,08275 0,07258 0,02717

#>

#> Obrót (n x k) = (21 x 21):

#> PC1 PC2 PC3 PC4 PC5

#> ID 0,008492 -0,007276 0,14499 0,174484 -0,82840

#> Mieszkańcy 0,205721 0,004321 0,54743 0,303663 0,06659

#> Gospodarstwa domowe 0,181071 0,008752 0,49902 0,470793 0,13119

#> AdultMeanAge -0,275210 0,192311 0,14601 -0,011834 0,12951

#> Biały -0,239842 0,112711 -0,25766 0,471189 -0,02500

#> Posiadane -0,289544 0,085502 0,26954 -0,179515 -0,11673

(Obcięte wyjście)

Kiedy wyświetlamy obiekt pca, możemy zobaczyć odchylenia standardowe dla każdej zmiennej, ale co ważniejsze, możemy zobaczyć wagi używane dla każdej zmiennej do utworzenia każdego głównego składnika. Jak widzimy, gdy spojrzymy na pełną wydajność naszego komputera, wśród najważniejszych wag (największych wartości bezwzględnych) mamy zmienne wieku i pochodzenia etnicznego, a także inne, takie jak posiadanie domu. Jeśli chcesz uzyskać wartość osi dla każdej obserwacji w nowym układzie współrzędnych złożonym z głównych składowych, po prostu musisz pomnożyć każdą obserwację w swoich danych (każdy wiersz) przez odpowiednie wagi z macierzy rotacji z obiektu pca (pca$rotation). Na przykład, aby wiedzieć, gdzie należy umieścić pierwszą obserwację w danych w odniesieniu do drugiego głównego składnika, możesz skorzystać z:

as.matrix(data_numerical_adjusted[1, ]) %*% pca$rotation [, 1]

Ogólnie rzecz biorąc, możesz zastosować operacje na macierzach, aby uzyskać współrzędne dla wszystkich obserwacji w danych w odniesieniu do wszystkich głównych składników w obiekcie 􀁑􀁄􀁂, używając następującego wiersza, który wykona mnożenie macierzy. Pamiętaj, że nie musisz tego robić samodzielnie, ponieważ R zrobi to automatycznie podczas analizy wyników.

as.matrix(data_numerical_adjusted) %*% pca$rotation

Kiedy patrzymy na podsumowanie pca, możemy zobaczyć odchylenia standardowe dla każdego głównego składnika, a także jego proporcję uchwyconej wariancji i jej akumulację. Informacje te są przydatne przy podejmowaniu decyzji, ile głównych składników powinniśmy zachować do końca analizy. W naszym przypadku okazuje się, że mając tylko pierwsze dwa główne składniki, uchwyciliśmy około 70 procent informacji zawartych w danych, co w naszym przypadku może być wystarczające. Liczbę 70% można uzyskać, dodając wartość Proportion of variance dla głównych składników, które chcemy rozważyć (w kolejności i zaczynając od PC1). W tym przypadku, jeśli dodamy Proportion of variance dla PC1 i PC2, otrzymamy 0,411$ + 0,280 = 0,691$, czyli prawie 70 proc. Zauważ, że możesz po prostu spojrzeć na Cumulative proportion do znalezienia tej liczby bez konieczności samodzielnego obliczania sumy, ponieważ gromadzi ona Proportion of variance stopniowo, zaczynając od PC1. Zastanów się przez chwilę, jak potężna jest ta technika: mając tylko dwie zmienne, jesteśmy w stanie uchwycić 70 procent informacji zawartych w pierwotnych 40 zmiennych:

summary(pca)

#> Znaczenie komponentów:

#> PC1 PC2 PC3 PC4 PC5 PC6 PC7

#> Odchylenie standardowe 2,939 2,426 1,2586 1,1330 1,0080 0,9411 0,7139

#> Odsetek wariancji 0,411 0,280 0,0754 0,0611 0,0484 0,0422 0,0243

#> Skumulowana proporcja 0,411 0,692 0,7670 0,8281 0,8765 0,9186 0,9429

#> PC8 PC9 PC10 PC11 PC12 PC13

#> Odchylenie standardowe 0,5761 0,5405 0,44767 0,37701 0,30166 0,212 11

#> Odsetek wariancji 0,0158 0,0139 0,00954 0,00677 0,00433 0,00214

#> Skumulowana proporcja 0,9587 0,9726 0,98217 0,98894 0,99327 0,99541

#> PC14 PC15 PC16 PC17 PC18 PC19

#> Odchylenie standardowe 0,17316 0,1376 0,11474 0,10843 0,09797 0,08275

#> Proporcja wariancji 0,00143 0,0009 0,00063 0,00056 0,00046 0,00033

#> Łączna proporcja 0,99684 0,9977 0,99837 0,99893 0,99939 0,99971

(Obcięte wyjście)

Na powyższym wykresie widzimy wariancje (w postaci kwadratów odchyleń standardowych) z wyników summary(pca). Widzimy, jak każdy kolejny składnik główny oddaje mniejszą wartość całkowitej wariancji:

plot(pca, type = „1” , main = „Principal Components’ Variaces”)

Wreszcie, poniższy wykres przedstawia wykres rozrzutu obserwacji (punktów) oddziału na płaszczyźnie utworzonej przez dwa główne składniki z naszej analizy; nazywa się to biplot.

 Ponieważ te dwa główne komponenty są tworzone jako liniowe kombinacje oryginalnych zmiennych, potrzebujemy pewnych wskazówek przy ich interpretacji. Dla ułatwienia strzałki wskazują kierunek przyporządkowania tej zmiennej do głównej osi składowej. Im dalej strzałka od środka, tym silniejszy wpływ na główne składniki. Dzięki temu biplotowi widzimy, że Proportion jest silnie powiązany z wychowankami, którzy głosowali za opuszczeniem UE, co jest oczywiste, ponieważ jest to konstrukcyjne. Jednak widzimy też inne interesujące relacje. Na przykład, poza efektami, które do tej pory stwierdziliśmy (wiek, wykształcenie i pochodzenie etniczne), osoby posiadające własne domy również są nieco powiązane z większą skłonnością do głosowania za opuszczeniem UE. Z drugiej strony nieznaną wcześniej relacją jest fakt, że im gęstsza jest populacja podopiecznego (pomyśl o miastach), tym bardziej prawdopodobne jest, że będą głosować za pozostaniem w UE:

library(ggbiplot)

biplot <- ggbiplot(pca, groups – data$Vote)

biplot <- biplot + scale_color_disrete(name = „ „)

biplot <- biplot + theme(legend.position = „top”, legend.discretio = „horizontal”)

print(biplot)

Tworzenie nowego zestawu danych na podstawie tego, czego się nauczyliśmy

Dowiedzieliśmy się dotychczas, że wiek, wykształcenie i pochodzenie etniczne są ważnymi czynnikami w zrozumieniu sposobu, w jaki ludzie głosowali w referendum w sprawie Brexitu. Młodsze osoby z wyższym wykształceniem kojarzone są z głosami za pozostaniem w UE. Starsi biali ludzie są kojarzeni z głosami za opuszczeniem UE. Możemy teraz wykorzystać tę wiedzę, aby stworzyć bardziej zwięzły zestaw danych, który obejmuje tę wiedzę. Najpierw dodajemy odpowiednie zmienne, a następnie usuwamy nieistotne zmienne. Nasze nowe istotne zmienne to dwie grupy wiekowe (dorośli poniżej i powyżej 45 lat), dwie grupy etniczne (biali i nie-biali) oraz dwie grupy edukacji (wysoki i niski poziom wykształcenia):

data$Age_18to44 <- (

data$Age_18to19 +

data$Age_20to24  +

data$Age_25to29 +

data$Age_30to44

)

data$Age_45plus <- (

data$Age_45to59 +

data$Age_60to64 +

data$Age_65to74 +

data$Age_75to84 +

data$Age_85to89 +

data$Age_90plus

)

data$NonWhite <- (

data$Black +

data$Asian +

data$Indian +

data$Pakistani

)

data$HighEducationLevel <- data$L4Quals_plus

data$LowEducationLevel1 <- data$oQuals

Teraz usuwamy stare zmienne, które były używane do tworzenia naszych nowo dodanych zmiennych. Aby to zrobić bez konieczności ręcznego określania pełnej listy, wykorzystując fakt, że wszystkie zawierają słowo „Age”, tworzymy wektor logiczny age_variables który zawiera wartość TRUE dla tych zmiennych, które zawierają słowo „Age” wewnątrz (w przeciwnym razie FALSE) i upewnij się, że zachowujemy nowo utworzone  zmienne Age_18to44  i Age_45plus. Pozostałe poziomy etniczności i wykształcenia usuwamy ręcznie:

column_names <- colnames(data)

new_variables <- !logical(length(columns_names))

new_variables <- setNames(new_variables, colimn_names)

age_variables <- sapply(column_names, function(x) grepl(„Age”, x))

new_variables [age_variables= <- FALSe

new_variables [[„AdultMeanAge”]] <- TRUE

new_variables [[„Age_18to44”]] <- TRUE

new_variables [[„Age_45plus”]] <- TRUE

new_variables [[„Black”]] <- FALSE

new_variables [[„Asian”]] <- FALSE

new_variables [[„India”]] <- FALSE

new_variables [[„Pakistani”]] <- FALSE

new_variables [[„NoQuals”]] <- FALSE

new_variables [[„L4Quals_plus”]] <- FALSE

new_variables [[„OwnedOutright”]] <- FALSE

new_variables [[„MultiDeprived”]] <- FALSE

Zapisujemy utworzony obiekt data_adjusted, wybierając nowe kolumny, tworzymy nowe zmienne numeryczne dla nowej struktury danych i zapisujemy go jako plik CSV:

data_adjusted <- data[, new_variables]

numerical_variables_adjusted <- sapply(data_adjusted, is.numeric)

write.csv(data_adjusted, file = „data_brexit_referendum_adjusted.csv”)

Zrozumienie interakcji z korelacjami

Korelacja jest miarą liniowej zależności między dwiema zmiennymi. Jego wartość waha się od -1, oznaczającego doskonałą relację odwrotną, do 1, reprezentującego doskonały związek bezpośredni. Podobnie jak utworzyliśmy macierz wykresów punktowych, teraz utworzymy macierz korelacji, a wynikowy wykres pokazano poniżej. Duże koła oznaczają wysoką korelację bezwzględną. Niebieskie kółka oznaczają korelację dodatnią, a kółka czerwone korelację ujemną. Aby utworzyć ten wykres, użyjemy funkcji corrplot() z pakietu corrplot i przekażemy jej dane korelacji obliczone przez funkcję cor() w R, i opcjonalnie niektóre parametry etykiet tekstowych (tl), takie jak kolor (color) i rozmiar (cex).

Spójrzmy teraz na następujący kod:

library(corrplot)

corrplot(corr = cor(data_numerical(, tl.col = „black”, tl.cex = 0.6)

Jeśli spojrzymy na relację między zmienną Proportion a innymi zmiennymi, zmienne w dużych niebieskich kółkach są z nią dodatnio skorelowane, co oznacza, że ​​im bardziej ta zmienna rośnie, tym większe jest prawdopodobieństwo, że zmienna Proportion również wzrosnie. Aby zobaczyć przykłady tego typu, spójrz na relacje między AdultMeanAge i NoQuals z Proportion. Jeśli wśród Proportion i innych zmiennych znajdziemy duże czerwone kółka, oznacza to, że tym bardziej zmienna rośnie, tym bardziej prawdopodobne jest, że Proportion spadnie. Aby zobaczyć przykłady tego typu, spójrz na relacje między Age_25to29, Age_30to44 i L4Quals+plus z Proportion:

Uzyskanie lepszego wyglądu dzięki szczegółowym wykresom punktowym

Teraz, gdy wiemy, jak uzyskać ogólny obraz wykresów rozrzutu, aby uzyskać ogólny sens relacji między zmiennymi, jak możemy uzyskać bardziej szczegółowe spojrzenie na każdy wykres punktowy? Cóż, cieszę się, że zapytałeś! Aby to osiągnąć, zrobimy to w dwóch krokach. Najpierw będziemy pracować nad stworzeniem pojedynczego, szczegółowego wykresu punktowego, z którego jesteśmy zadowoleni. Po drugie, opracujemy prosty algorytm, który przejdzie przez wszystkie kombinacje zmiennych i utworzy odpowiedni wykres dla każdej z nich:

Wykres pokazany powyżej przedstawia nasz prototypowy wykres punktowy. Ma kombinację zmiennych w osiach x i y, NoQuals i AdultMeanAge w naszym przypadku przypisuje kolor zgodnie z odpowiednim Proportion i umieszcza na górze linię odpowiadającą regresji liniowej, aby uzyskać ogólny obraz relacji między zmiennymi na osiach. Porównaj ten wykres z lewym wykresem punktowym z poprzedniej pary wykresów punktowych. To ta sama fabuła, ale ta jest bardziej szczegółowa i zawiera więcej informacji. Ta fabuła wydaje się na razie wystarczająco dobra.

plot ,- ggplot(data, aes(x = NoQuals, y = AdultMeanAge, color = Proportion))

plot <- plot + stat_smooth(method = „lm”, col = „darkgrey”,  se= FALSE)

plot <- plot + scale_color_viridis()

plot <- plot  + geom_point()

print(plot)

Teraz musimy opracować algorytm, który weźmie wszystkie kombinacje zmiennych i utworzy odpowiednie wykresy. Przedstawiamy pełny algorytm i wyjaśniamy część po części. Jak widać, zaczynamy definiować funkcję create_graphs_iteratively, która otrzymuje dwa parametry: data i plot_function. Algorytm pobierze nazwy zmiennych dla danych i zapisze je w zmiennych vars. Następnie usunie Proportio z takich zmiennych, ponieważ będą one używane do tworzenia kombinacji dla osi, a Proportio nigdy nie zostanie użyte na osi; będzie używany wyłącznie do kolorów. Teraz, jeśli wyobrazimy sobie wszystkie kombinacje zmiennych w macierzy, takiej jak ta dla wykresu punktowego macierzy pokazanego wcześniej, to musimy przejść przez górny lub dolny trójkąt, aby uzyskać wszystkie możliwe kombinacje (w rzeczywistości górny i dolny trójkąt z macierze wykresów punktowych są symetryczne, ponieważ przekazują te same informacje). Aby przejść przez te trójkąty, możemy użyć znanego wzoru, który wykorzystuje dwie pętle for, każda dla jednej osi, i gdzie wewnętrzna pętla musi zaczynać się tylko w miejscu zewnętrznej pętli (to tworzy trójkąt). -1 i +1 są po to, aby upewnić się, że zaczynamy i kończymy w odpowiednich miejscach w każdej pętli bez otrzymania błędu dotyczącego granic tablicy. Wewnątrz pętli wewnętrznej utworzymy nazwę wykresu jako kombinację nazw zmiennych i połączymy je za pomocą funkcji paste(), a także utworzymy wykres za pomocą plot_fuction wyślemy jako parametr (więcej na ten temat w dalszej części). Funkcje pg() i dev.off() służą do zapisywania wykresów na dysku twardym komputera. Pomyśl o funkcji png() jako miejscu, w którym R zaczyna szukać wykresu, a dev.off() jako miejscu, w którym zatrzymuje proces zapisywania. Zapraszam do zajrzenia do ich dokumentacji lub poczytać więcej o urządzeniach w R.

create_plots+iteratively <- function(data, plot_function) {

vars <- colnames(data)

vars <- vars(!which(vars ==”Proportion”))

for (i in 1 : (length(vars) – 1 )) {

for (j in (i + 1) : length (vars) ) {

save_to <- paste(vars[i], „ _ „ , vars[j], „.png”, sep = „ „)

plot_function(data, vars[i[, vars[j], save_to)

}

}

}

Prawie skończyliśmy; musimy tylko opakować kod, którego użyliśmy do przekształcenia naszego prototypu wykresu w funkcję i wszystko będzie gotowe. Jak widać, wyodrębniliśmy parametry x, y i color dla funkcji aes() jako zmienne, które są wysyłane jako parametry do funkcji (nazywa się to argumentami parametryzującymi), a zmieniliśmy funkcję aes()z funkcją aes_string(), która może otrzymywać zmienne z łańcuchami parametrów. Dodaliśmy również opcję wysyłania var_color jako FALSE, aby uniknąć używania kolorowej wersji wykresu. Wszystko inne jest takie samo:

prototype_scatter_plot <- function(data, var_x, var_y, var_color = „Proportion”, save_to = „ „ ) {

if (is.na(as.logical(var_color))) {

plot <- ggplot(data, aes_string(x = var_x, y = var_y, color = var_color))

} else{

plot <- ggplot(data, aes_string(x = var_x, y = var_y))

}

plot <- plot + stat_smooth (method = „lm”, col = „darkgrey”, se = FALSE)

plot <- plot + scale_color_viridis()

plot <- plot + geom_point()

if (not_empty(save_to)) png(save_to)

print(plot)

if (not_empty(save_to)) dev.off

}

Ponieważ będziemy sprawdzać w różnych miejscach, czy łańcuch save_to jest pusty, nazywamy czek i zawijamy go w funkcję not_empty(). Teraz trochę łatwiej jest czytać nasz kod.

not_empty <- function(file) {

return(file != „”)

}

Dzięki tej funkcji protoype_scatter_plot() możemy dość łatwo odtworzyć prawidłowe wykresy rozrzutu pokazane wcześniej, a także dowolną inną kombinację zmiennych . Wydaje się to dość potężne, prawda?

Spójrzmy na następujący kod:

protoype_scatter_plot(data, „L4Quals_plus”, „AdultMeanAge”

Teraz, gdy wykonaliśmy ciężką pracę, możemy dość łatwo stworzyć wszystkie możliwe kombinacje. Musimy tylko wywołać funkcję create_plots_iteratively z naszymi danymi i funkcja prototype_scatter_plot(). Używanie funkcji jako parametrów dla innych funkcji jest znane jako wzorzec strategii. Nazwa wzięła się stąd, że możemy łatwo zmienić naszą strategię tworzenia wykresów dla każdego innego, który chcemy, a który otrzymuje te same parametry (data, var_x i var_y), aby tworzyć wykresy, bez konieczności zmiany naszego algorytmu, aby przechodzić przez kombinacje zmiennych. Ten rodzaj elastyczności jest bardzo potężny:

create_plots_iteratively(data, prototype_scatter_plot)

Spowoduje to utworzenie dla nas wszystkich działek i zapisanie ich na naszym dysku twardym. Całkiem fajnie, co? Teraz możemy spojrzeć na każdy z nich niezależnie i zrobić z nimi wszystko, co potrzebujemy, ponieważ mamy je już jako pliki PNG.

Użycie macierzowych wykresów punktowych dla szybkiego przeglądu

Co się stanie, jeśli chcemy zwizualizować wiele wykresów rozrzutu na jednym wykresie, aby szybko uzyskać sens w danych? W takim przypadku potrzebujemy macierzowych wykresów punktowych. Mamy różne opcje pakietu do tworzenia takich macierzowych wykresów punktowych (np. pakiet car). Jednak dla uproszczenia użyjemy funkcji wbudowanej zamiast pakietu zewnętrznego. Spoglądając na poniższy wykres, możemy uzyskać ogólny obraz interakcji między zmiennymi. Celem tego typu wizualizacji nie jest przedstawienie szczegółów, ale przedstawienie ogólnego przeglądu. Aby przeczytać ten wykres, musimy spojrzeć na dowolny interesujący wykres punktowy w macierzy i przesuwać się zarówno w poziomie, jak iw pionie, aż znajdziemy nazwę związaną z jego osią. Na przykład, jeśli spojrzysz na wykres bezpośrednio po prawej stronie NoQuals i jednocześnie bezpośrednio nad L4Quals_plus, to patrzysz na relację między tymi dwiema zmiennymi (NoQuals na osi y, L4Quals_plus,  na osi x) i okazuje się, że jest to relacja odwrotna; im wyższy odsetek osób na oddziale o wysokim poziomie wykształcenia, tym niższy odsetek osób z niskim wykształceniem. Inną oczywistą zależnością jest to, że im wyższy poziom wykształcenia (L4Quals_plus,), tym wyższy zawód (HigherOccup).

Ze względu na ograniczenia przestrzenne nie byliśmy w stanie pokazać wszystkich relacji zmiennych, ponieważ wykresy rozrzutu byłyby zbyt małe, aby nadać im sens. Zachęcamy jednak czytelnika do dodawania kolejnych zmiennych do macierzy. Istnieją pewne nieoczywiste relacje. Znalezienie ich pozostawiamy jako ćwiczenie dla czytelnika:

desired_variablees <- c)

„AdultMeanAge”,

„White”,

„Owned”,

„NoQuals”

„L4Quals_plus”,

„Unemp:,

„HigherOccup”,

„Deprived”;

„Proportion”

)

pairs(data[,desired_variables])