Statystyka matematyczna (kierunek informatyka i ekonometria)

Praca z danymi i testowanie danych

¹ Zakład Bioinformatyki, Instytut Informatyki, Uniwersytet w Białymstoku

^✉ Correspondence: Jarosław Kotowicz <j.kotowicz@uwb.edu.pl>

1 Grupa 2 (6 maja 2020r.)

1.1 Praca z danymi i testowanie hipotez

1.1.1 Czyszczenie danych globalnych i podczytanie danych

rm(list = ls())

library(readr)
dane <- read_csv("http://www.jaredlander.com/data/acsNew.csv")

Parsed with column specification:
cols(
  Acres = col_character(),
  FamilyIncome = col_double(),
  FamilyType = col_character(),
  NumBedrooms = col_double(),
  NumChildren = col_double(),
  NumPeople = col_double(),
  NumRooms = col_double(),
  NumUnits = col_character(),
  NumVehicles = col_double(),
  NumWorkers = col_double(),
  OwnRent = col_character(),
  YearBuilt = col_character(),
  HouseCosts = col_double(),
  ElectricBill = col_double(),
  FoodStamp = col_character(),
  HeatingFuel = col_character(),
  Insurance = col_double(),
  Language = col_character(),
  Income = col_character()
)
Ostrzeżenie: zamykanie nieużywanego połączenia 3 (http://www.jaredlander.com/data/acsNew.csv)

1.1.2 Praca z danymi (tworzenie nowych zmiennych)

library(tidyverse)

Registered S3 method overwritten by 'dplyr':
  method           from
  print.rowwise_df     
Registered S3 methods overwritten by 'dbplyr':
  method         from
  print.tbl_lazy     
  print.tbl_sql      
-- Attaching packages --------------------------------------- tidyverse 1.3.0 --
<U+221A> ggplot2 3.3.0     <U+221A> dplyr   0.8.5
<U+221A> tibble  3.0.1     <U+221A> stringr 1.4.0
<U+221A> tidyr   1.0.2     <U+221A> forcats 0.5.0
<U+221A> purrr   0.3.4     
-- Conflicts ------------------------------------------ tidyverse_conflicts() --
x dplyr::filter() masks stats::filter()
x dplyr::lag()    masks stats::lag()

dane <- dane %>%
  mutate(FamilyIncomeSd = (FamilyIncome - mean(FamilyIncome))/sd(FamilyIncome))

1.1.3 Testowanie hipotez nieparametrycznych

1.1.3.1 Testowanie zgodności z rozkładem normalnym

library(nortest)

Jak nie należy wywoływać (brak nazw zmiennych w środowisku)

ad.test(FamilyIncomesd)

Błąd w poleceniu 'sort(x[complete.cases(x)])':
  nie znaleziono obiektu 'FamilyIncomesd'

Prawidłowy sposób wywołania

ad.test(dane$FamilyIncomeSd)

    Anderson-Darling normality test

data:  dane$FamilyIncomeSd
A = 132.37, p-value < 2.2e-16

Interpretacja wyniku!

Można dołączyć nazwy zmiennych z ramki danych do zmiennych globalnych

attach(dane)

I wywołać test

ad.test(FamilyIncomeSd)

    Anderson-Darling normality test

data:  FamilyIncomeSd
A = 132.37, p-value < 2.2e-16

Interpretacja wyniku!

Usunięcie nazwy zmiennych ze środowiska

detach(dane)

Ponownie nie działa

ad.test(FamilyIncomesd)

Błąd w poleceniu 'sort(x[complete.cases(x)])':
  nie znaleziono obiektu 'FamilyIncomesd'

Test Kołomogorowa-Smirnowa

ks.test(dane$FamilyIncomeSd, "pnorm")

warto㤼㹣ci powt昼㸳rzone nie powinny by攼㸶 obecne w te㤼㹣cie Kolmogorowa-Smirnowa

    One-sample Kolmogorov-Smirnov test

data:  dane$FamilyIncomeSd
D = 0.15075, p-value < 2.2e-16
alternative hypothesis: two-sided

Interpretacja wyniku!

mean(dane$FamilyIncomeSd)

[1] 6.661647e-17

sd(dane$FamilyIncomeSd)

[1] 1

Interpretacja wyniku!

1.1.3.2 Generowanie ciągu liczb psudolosowych z rozkładu jednostajnego

set.seed(1115)
a <- runif(1000, -10, 5)

1.1.3.3 Testowanie zgodności z rozkładem jednostajnym

ks.test(a, "punif")

    One-sample Kolmogorov-Smirnov test

data:  a
D = 0.66, p-value < 2.2e-16
alternative hypothesis: two-sided

Interpretacja wyniku!

ks.test(a, "punif", -10, 5)

    One-sample Kolmogorov-Smirnov test

data:  a
D = 0.024471, p-value = 0.5872
alternative hypothesis: two-sided

Interpretacja wyniku!

ks.test(a, "punif", -12, 6)

    One-sample Kolmogorov-Smirnov test

data:  a
D = 0.11141, p-value = 3.313e-11
alternative hypothesis: two-sided

Interpretacja wyniku!

2 Grupa 1 (13 maja 2020r.)

2.1 Praca z danymi i testowanie hipotez

2.1.1 Czyszczenie danych globalnych i podczytanie danych

rm(list=ls())

# library(nortest)
# library(tidyverse)

dane <- read_csv("http://www.jaredlander.com/data/acsNew.csv")

Parsed with column specification:
cols(
  Acres = col_character(),
  FamilyIncome = col_double(),
  FamilyType = col_character(),
  NumBedrooms = col_double(),
  NumChildren = col_double(),
  NumPeople = col_double(),
  NumRooms = col_double(),
  NumUnits = col_character(),
  NumVehicles = col_double(),
  NumWorkers = col_double(),
  OwnRent = col_character(),
  YearBuilt = col_character(),
  HouseCosts = col_double(),
  ElectricBill = col_double(),
  FoodStamp = col_character(),
  HeatingFuel = col_character(),
  Insurance = col_double(),
  Language = col_character(),
  Income = col_character()
)

dane <- dane %>%
  mutate_if(is.character,list(as.factor))

dane <- dane %>%
  mutate(FamilyIncomeSd = (FamilyIncome - mean(FamilyIncome))/sd(FamilyIncome))

2.1.1.1 Dołącznie zmiennych (nazw z ramki danych) do środowiska

attach(dane)

2.1.2 Testy zgodności z rozkładem normalnym

sf.test(FamilyIncomeSd)

    Shapiro-Francia normality test

data:  FamilyIncomeSd
W = 0.73234, p-value < 2.2e-16

Interpretacja wyniku!

detach(dane)

2.1.3 Generowanie ciągu liczb pseudolosowych

set.seed(1050)
x.norm01 <- rnorm(1000)

2.1.4 Ponownie testy zgodności z rozkładem normalnym

cvm.test(x.norm01)

    Cramer-von Mises normality test

data:  x.norm01
W = 0.039185, p-value = 0.6967

Interpretacja wyniku!

2.1.4.1 Test zgodności Kołmogorowa-Smirnowa (zgodności z rozkładem wzorocwym i zgodności dóWch rozkładów)

Jak nie robić! (podczytuje wartości domyślne)

ks.test(x.norm01, "pnorm")

    One-sample Kolmogorov-Smirnov test

data:  x.norm01
D = 0.037062, p-value = 0.1282
alternative hypothesis: two-sided

Interpretacja wyniku!

set.seed(1050)
x.norm32 <- rnorm(1000, 3 ,2)

Interpretacja wyniku!

ks.test(x.norm32, "pnorm",0 ,1)

    One-sample Kolmogorov-Smirnov test

data:  x.norm32
D = 0.74095, p-value < 2.2e-16
alternative hypothesis: two-sided

Interpretacja wyniku!

ks.test(x.norm32, "pnorm", 3, 2)

    One-sample Kolmogorov-Smirnov test

data:  x.norm32
D = 0.037062, p-value = 0.1282
alternative hypothesis: two-sided

Interpretacja wyniku!

set.seed(1050)
x.unif <- runif(1000, -1, 1)

ks.test(x.unif, "punif")

    One-sample Kolmogorov-Smirnov test

data:  x.unif
D = 0.50657, p-value < 2.2e-16
alternative hypothesis: two-sided

Interpretacja wyniku!

ks.test(x.unif, "punif",-1, 1)

    One-sample Kolmogorov-Smirnov test

data:  x.unif
D = 0.016202, p-value = 0.9555
alternative hypothesis: two-sided

Interpretacja wyniku!

ks.test(x.unif[1:100],x.unif[901:1000])

    Two-sample Kolmogorov-Smirnov test

data:  x.unif[1:100] and x.unif[901:1000]
D = 0.12, p-value = 0.4676
alternative hypothesis: two-sided

Interpretacja wyniku!

ks.test(x.unif[1:100],x.unif[901:1000],alternative = "less")

    Two-sample Kolmogorov-Smirnov test

data:  x.unif[1:100] and x.unif[901:1000]
D^- = 0.07, p-value = 0.6126
alternative hypothesis: the CDF of x lies below that of y

Interpretacja wyniku!

ks.test(x.unif[1:100],x.unif[901:1000],alternative = "greater")

    Two-sample Kolmogorov-Smirnov test

data:  x.unif[1:100] and x.unif[901:1000]
D^+ = 0.12, p-value = 0.2369
alternative hypothesis: the CDF of x lies above that of y

2.1.5 Testy parametryczne równości wariancji

var.test(x.norm01, x.norm32)

    F test to compare two variances

data:  x.norm01 and x.norm32
F = 0.25, num df = 999, denom df = 999, p-value < 2.2e-16
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
 0.2208247 0.2830300
sample estimates:
ratio of variances 
              0.25 

Interpretacja wyniku!

ansari.test(x.norm01, x.norm32)

    Ansari-Bradley test

data:  x.norm01 and x.norm32
AB = 537652, p-value = 8.709e-09
alternative hypothesis: true ratio of scales is not equal to 1

Interpretacja wyniku!

mood.test(x.norm01[1:100], x.norm01[901:1000])

    Mood two-sample test of scale

data:  x.norm01[1:100] and x.norm01[901:1000]
Z = 0.98507, p-value = 0.3246
alternative hypothesis: two.sided

Interpretacja wyniku!

2.2 Regresja

2.2.1 Budowa modelu

model_01 <- lm(FamilyIncome ~ NumPeople, dane)

2.2.2 Podsumowanie modelu

summary(model_01)

Call:
lm(formula = FamilyIncome ~ NumPeople, data = dane)

Residuals:
    Min      1Q  Median      3Q     Max 
-150196  -55913  -22018   25329  905682 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)    88377       5336  16.563  < 2e-16 ***
NumPeople       6647       1450   4.583 4.84e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 96970 on 2271 degrees of freedom
Multiple R-squared:  0.009163,  Adjusted R-squared:  0.008727 
F-statistic:    21 on 1 and 2271 DF,  p-value: 4.836e-06

Interpretacja wyniku!

2.2.3 Wykresy diagnostyczne dla modelu

plot(model_01)

2.2.4 Drugi model

model_02 <- lm(FamilyIncomeSd ~ NumPeople, dane)
summary(model_02)

Call:
lm(formula = FamilyIncomeSd ~ NumPeople, data = dane)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.5421 -0.5741 -0.2261  0.2601  9.2986 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.23209    0.05478  -4.237 2.36e-05 ***
NumPeople    0.06825    0.01489   4.583 4.84e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9956 on 2271 degrees of freedom
Multiple R-squared:  0.009163,  Adjusted R-squared:  0.008727 
F-statistic:    21 on 1 and 2271 DF,  p-value: 4.836e-06

Interpretacja wyniku!

detach(package:tidyverse)
detach(package:ggplot2)
detach(package:tibble)
detach(package:tidyr)
detach(package:readr)
detach(package:purrr)
detach(package:dplyr)
detach(package:stringr)
detach(package:forcats)

detach(package:nortest)