R za novu generaciju IT stručnjaka

U današnjem ćemo tekstu prezentirati primjer jednostavne analize podataka unutar sustava R. Podaci koje ćemo ovom prilikom koristiti dolaze zajedno s instalacijom sustava i sadržavaju skup podataka o cijenama prodanih kuća u gradu Amesu, Iowa u razdoblju od 2006. do 2010.

Dodatno, skup podataka moguće je pronaći i na drugim web stranicama poput http://lib.stat.cmu.edu/datasets/boston. Skup podataka sadrži brojne varijable o kvaliteti i količini fizičkih atributa stambenih kuća u Iowi prodanih u razdoblju od 2006. do 2010. godine. To su najčešće varijable koje kupac želi znati o nekretnini (kvadratura, broj spavaćih soba i kupaonica, veličina parcele itd.).

Kako ovaj skup podataka dolazi instalacijom sustava R, dovoljno je u konzolu upisati ime objekta (skupa) kako bi se ispisao na konzoli. Također, sve funkcije sustava možemo odmah primijeniti na podatke.

Pogledajmo prva dva retka skupa podataka.

ames[1:2,]

Order PID MS.SubClass MS.Zoning Lot.Frontage Lot.Area Street Alley

1 1 526301100 20 RL 141 31770 Pave

2 2 526350040 20 RH 80 11622 Pave

Lot.Shape Land.Contour Utilities Lot.Config Land.Slope Neighborhood

1 IR1 Lvl AllPub Corner Gtl NAmes

2 Reg Lvl AllPub Inside Gtl NAmes

Condition.1 Condition.2 Bldg.Type House.Style Overall.Qual Overall.Cond

1 Norm Norm 1Fam 1Story 6 5

2 Feedr Norm 1Fam 1Story 5 6

Year.Built Year.Remod.Add Roof.Style Roof.Matl Exterior.1st Exterior.2nd

1 1960 1960 Hip CompShg BrkFace Plywood

2 1961 1961 Gable CompShg VinylSd VinylSd

Mas.Vnr.Type Mas.Vnr.Area Exter.Qual Exter.Cond Foundation Bsmt.Qual

1 Stone 112 TA TA CBlock TA

2 None 0 TA TA CBlock TA

Bsmt.Cond Bsmt.Exposure BsmtFin.Type.1 BsmtFin.SF.1 BsmtFin.Type.2

1 Gd Gd BLQ 639 Unf

2 TA No Rec 468 LwQ

BsmtFin.SF.2 Bsmt.Unf.SF Total.Bsmt.SF Heating Heating.QC Central.Air

1 0 441 1080 GasA Fa Y

2 144 270 882 GasA TA Y

Electrical X1st.Flr.SF X2nd.Flr.SF Low.Qual.Fin.SF Gr.Liv.Area

1 SBrkr 1656 0 0 1656

2 SBrkr 896 0 0 896

Bsmt.Full.Bath Bsmt.Half.Bath Full.Bath Half.Bath Bedroom.AbvGr

1 1 0 1 0 3

2 0 0 1 0 2

Kitchen.AbvGr Kitchen.Qual TotRms.AbvGrd Functional Fireplaces

1 1 TA 7 Typ 2

2 1 TA 5 Typ 0

Fireplace.Qu Garage.Type Garage.Yr.Blt Garage.Finish Garage.Cars

1 Gd Attchd 1960 Fin 2

2 Attchd 1961 Unf 1

Garage.Area Garage.Qual Garage.Cond Paved.Drive Wood.Deck.SF

1 528 TA TA P 210

2 730 TA TA Y 140

Open.Porch.SF Enclosed.Porch X3Ssn.Porch Screen.Porch Pool.Area Pool.QC

1 62 0 0 0 0

2 0 0 0 120 0

Fence Misc.Feature Misc.Val Mo.Sold Yr.Sold Sale.Type Sale.Condition

1 0 5 2010 WD Normal

2 MnPrv 0 6 2010 WD Normal

SalePrice

1 215000

2 105000

Tijekom analize postupno ćemo se upoznavati s podacima putem jednostavnih pitanja kao što je: “Koja je dimenzija podataka?”

dim(ames)

[1] 2930 82

Koliko varijabli se nalazi u skupu podataka?

dim(ames)[2] #ili ncol(ames)

[1] 82

Koliko opservacija (promatranja) sadrži skup podataka?

dim(ames)[1] # ili nrow(ames)

[1] 2930

Koja je klasa objekta?

class(ames)

[1] “data.frame”

Imena varijabli u skupu podataka.

names(ames)

[1] “Order” “PID” “MS.SubClass”

[4] “MS.Zoning” “Lot.Frontage” “Lot.Area”

[7] “Street” “Alley” “Lot.Shape”

[10] “Land.Contour” “Utilities” “Lot.Config”

[13] “Land.Slope” “Neighborhood” “Condition.1”

Sumarna statistika po varijablama – ispis za prvih 6 varijabli.

summary(ames)[,1:6]

Order PID MS.SubClass MS.Zoning

Min. : 1.0 Min. :5.263e+08 Min. : 20.00 A (agr): 2

1st Qu.: 733.2 1st Qu.:5.285e+08 1st Qu.: 20.00 C (all): 25

Median :1465.5 Median :5.355e+08 Median : 50.00 FV : 139

Mean :1465.5 Mean :7.145e+08 Mean : 57.39 I (all): 2

3rd Qu.:2197.8 3rd Qu.:9.072e+08 3rd Qu.: 70.00 RH : 27

Max. :2930.0 Max. :1.007e+09 Max. :190.00 RL :2273

RM : 462

Lot.Frontage Lot.Area

Min. : 21.00 Min. : 1300

1st Qu.: 58.00 1st Qu.: 7440

Median : 68.00 Median : 9436

Mean : 69.22 Mean : 10148

3rd Qu.: 80.00 3rd Qu.: 11555

Max. :313.00 Max. :215245

NA’s :490

Koliki je interkvartilni raspon varijable SalePrice – iznosa za koji je prodana nekretnina?

IQR(ames$SalePrice)

[1] 84000

Kolike su minimalne i maksimalne cijene nekretnina prodanih u gradu Amesu od 2006. do 2010. godine?

min(ames$SalePrice)

[1] 12789

max(ames$SalePrice)

[1] 755000

Koliko različitih vrijednosti nalazimo u varijabli Garage.Cond? Za koliko opservacija nemamo podatak u varijabli?

summary(ames$Garage.Cond)

Ex Fa Gd Po TA NA’s

1 3 74 15 14 2665 158

Iz ispisa sumarne statistike varijable vidljivo je da postoje 4 razine vrijednosti varijable te da za 158 opservacija u ovoj varijabli nema zapisa.

Koliki je raspon cijena nekretnina u ulicama koje imaju popločenje (Street == Pave)? Ulice su ili popločene ili pošljunčane (Pave/Grvl).

range(ames[ames$Street==”Pave”,]$SalePrice)

[1] 12789 755000

Nacrtat ćemo stem-and-leaf dijagram varijable SalePrice.

##

## The decimal point is 5 digit(s) to the right of the |

##

## 0 | 113444444

## 0 | 55555556666666666666666666666666777777777777777778888888888888888888+109

## 1 | 00000000000000000000000000000000000000000000000000000000000000000000+846

## 1 | 55555555555555555555555555555555555555555555555555555555555555555555+805

## 2 | 00000000000000000000000000000000000000000000000000000000000000000000+355

## 2 | 55555555555555555555555555555555555555555555555555555555555666666666+160

## 3 | 00000000000000000000000111111111111111111122222222222222222222222222+36

## 3 | 55555555555566666666666777777777777788888888888888889999999999

## 4 | 000000000111111112222222223334444

## 4 | 55556666777788899

## 5 | 00044

## 5 | 5566889

## 6 | 1123

## 6 |

## 7 |

## 7 | 56

Podatke ćemo prikazati i odgovarajućim grafikonima.

Cijena nekretnine – prikaz u 30 perioda podjele

Cijena nekretnine – prikaz prema godini prodaje

Koliko je kamina (Fireplaces) imala najskuplje prodana nekretnina u 2008. godini?

prodane_08

prodane_08[which.max(prodane_08$SalePrice),]$Fireplaces

[1] 2

Izračunat ćemo kovarijancu varijabli cijene kuće SalePrice i ukupne bruto površine Gr.Liv.Area. Kovarijanca je mjera linearne povezanosti dviju varijabli.

cov(ames$Gr.Liv.Area , ames$SalePrice)

[1] 28542200

Kako je kovarijanca mjera koja ovisi o mjernim jedinicama, ne govori nam mnogo tako dugo dok je ne standardiziramo, tj. izrazimo u obliku Pearson korelacijskog koeficijenta R.

cor( ames$SalePrice, ames$Gr.Liv.Area, method=”pearson”)

[1] 0.7067799

Približno 50 posto varijabilnosti cijene nekretnine (0.7 x 0.7) može se objasniti njezinom bruto površinom. Kako bismo jednostavnije odgovorili na ovakva pitanja, nacrtat ćemo dijagram raspršenja, engl. scatterplot, za varijable SalePrice (y os) i Gr.Liv.Area (x os). Izrada linearnog modela

Budući da smo zaključili da je veza naše varijable od interesa (cijena nekretnine) u linearnoj vezi s varijablom bruto površine, želimo ovu vezu i matematički opisati na cjelokupnom rasponu mogućih vrijednosti.

Kako ne bismo stalno pisali ime objekta (ames), napravit ćemo attach objekta i od sada pisati samo imena varijabli te izračunati jednadžbu pravca – linearni model.

attach(ames)

fit

summary(fit)

Call:

lm(formula = SalePrice ~ Gr.Liv.Area)

Residuals:

Min 1Q Median 3Q Max

-483467 -30219 -1966 22728 334323

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 13289.634 3269.703 4.064 4.94e-05 ***

Gr.Liv.Area 111.694 2.066 54.061

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1

Residual standard error: 56520 on 2928 degrees of freedom

Multiple R-squared: 0.4995, Adjusted R-squared: 0.4994

F-statistic: 2923 on 1 and 2928 DF, p-value:

Iz rezultata linearnog modela vidljivo je da je on statistički značajan ( p-value:

Interval pouzdanosti modela – predefinirano za 95 % interval ako ne tražimo drugačije.

confint(fit) #predefinirano za 0.95

2.5 % 97.5 %

(Intercept) 6878.4845 19700.7842

Gr.Liv.Area 107.6429 115.7451

Ako želimo procijeniti vrijednost neke nekretnine poznate bruto površine, trebali bismo izraditi, temeljem naših podataka i izmjerenoj kovarijanci, optimalan matematički model (pravac) kojim ćemo moći dati procjenu vrijednosti nekretnine za sve vrijednosti površina unutar kojih smo imali mjerenja u skupu podataka. Pravac se kroz skup točaka optimalno provlači metodom najmanjih kvadrata.

Ako želimo dobiti parametre pravca (linearnog modela), pogledamo objekt fit u koji smo spremili rezultat linearnog modeliranja, rezultat funkcije lm. Kako glasi jednadžba pravca kojim bruto površinom kuće objašnjavamo njezinu cijenu?

fit

Call:

lm(formula = SalePrice ~ Gr.Liv.Area)

Coefficients:

(Intercept) Gr.Liv.Area

13289.6 111.7

Vizualizirat ćemo linearni model (pravac) koji je optimalno provučen kroz točke te ga možemo koristiti za procjenu cijene nekretnine za bilo koju kvadraturu kuće. Procijenit ćemo cijenu nekretnine za 2000 kvadratnih feeta.

new

predict(fit, new,interval=”prediction”)

fit lwr upr

1 236677.6 125809 347546.2

Slično možemo napraviti i za neku drugu varijablu od interesa ili kombinacijom više varijabli. Primjerice, koliki postotak cijene kuće možemo objasniti površinom (Gr.Liv.Area) i brojem spavaćih soba (Bedroom.AbvGr).

fit2

summary(fit2)

Call:

lm(formula = SalePrice ~ Gr.Liv.Area + Bedroom.AbvGr)

Residuals:

Min 1Q Median 3Q Max

-581397 -27840 -862 24526 329875

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 59496.236 3741.249 15.90

Gr.Liv.Area 136.361 2.247 60.69

Bedroom.AbvGr -29149.110 1372.135 -21.24

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1

Residual standard error: 52620 on 2927 degrees of freedom

Multiple R-squared: 0.5664, Adjusted R-squared: 0.5661

F-statistic: 1912 on 2 and 2927 DF, p-value:

Stavljanjem dvije varijable kao objašnjavajuće, uspjeli smo objasniti 56,64 % varijabiliteta cijene nekretnina.

Zanima nas je li drugi model statistički bolji od modela samo s jednom varijablom. Provjerit ćemo omjere varijabiliteta (napraviti analizu varijance – ANOVA) koje objašnjavaju jedan i drugi model. Prisjetimo se da smo model samo s jednom varijablom (bruto površina) zvali fit, dok smo model s 2 varijable zvali fit2.

anova(fit, fit2, test=”F”)

Analysis of Variance Table

Model 1: SalePrice ~ Gr.Liv.Area

Model 2: SalePrice ~ Gr.Liv.Area + Bedroom.AbvGr

Res.Df RSS Df Sum of Sq F Pr(>F)

1 2928 9.3549e+12

2 2927 8.1052e+12 1 1.2497e+12 451.29

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1

Razlika u objašnjenom varijabilitetu dva modela je statistički veoma značajna – Pr(>F) 2.2e-16.

Dodavanjem većeg broja varijabli u model vrlo vjerojatno ćemo povećavati i stupanj objašnjenog varijabiliteta varijable od interesa, ali u nekom trenutku model će postati teško objašnjiv i kompleksan, ali i doveli do tzv. overfittinga. Zato je veoma važno od velikog skupa dostupnih varijabli (u našem slučaju 82) na optimalan način odabrati onaj podskup varijabli kojim ćemo dobiti što jednostavniji model, a koji pak objašnjava u najvećoj mjeri našu varijablu od interesa.

Metode kojima se radi selekcija podskupa varijabli kompleksan je dio statistike i obrade podataka. Radi demonstracije, napravit ćemo tzv. “stepwise” selekciju varijabli. U ovu proceduru trebali bismo staviti samo varijable koje imaju linearnu povezanost s varijablom od interesa, ali u ovoj demonstraciji stavljamo skup od 15 varijabli (redom od pozicije 10 do 25). Generalno, stepwise selekcijske metode mogu biti engl. forward, backward ili stepwise, a razlikuje ih način na koji stavljaju, odnosno izbacuju varijable koje ne doprinose značajno općoj kvaliteti modela i u ovome tekstu ih nećemo detaljno objašnjavati. Statistika koju u pozadini ove metode koriste je Akaike Information Criterion (AIC) – statistička mjera koja traži model s najmanjom Root Mean Squared Error (RMSE), ali pritom penalizira ukupan broj varijabli uključenih u model.

Postoje i druge metode selekcije podskupa varijabli u modelu koje koriste druge statistike kao indikatore optimalnog modela.

# model sa svim varijablama

full

# Stepwise regression model – u oba smjera

step

Start: AIC=61695.59

SalePrice ~ Land.Contour + Utilities + Lot.Config + Land.Slope +

Neighborhood + Condition.1 + Condition.2 + Bldg.Type + House.Style +

Overall.Qual + Overall.Cond + Year.Built + Year.Remod.Add +

Roof.Style + Roof.Matl + Exterior.1st

Df Sum of Sq RSS AIC

– Utilities 2 1.7908e+09 3.8313e+12 61693

3.8295e+12 61696

– Overall.Cond 1 4.6711e+09 3.8341e+12 61697

– Condition.2 7 3.1903e+10 3.8614e+12 61706

– Land.Contour 3 2.2329e+10 3.8518e+12 61707

– Condition.1 8 4.4214e+10 3.8737e+12 61713

– Lot.Config 4 4.1502e+10 3.8710e+12 61719

– Land.Slope 2 4.3752e+10 3.8732e+12 61725

– Year.Remod.Add 1 6.5731e+10 3.8952e+12 61743

– Year.Built 1 6.9731e+10 3.8992e+12 61746

– House.Style 7 8.8164e+10 3.9176e+12 61748

– Roof.Matl 7 1.0584e+11 3.9353e+12 61761

– Exterior.1st 15 1.3284e+11 3.9623e+12 61766

– Roof.Style 5 1.0921e+11 3.9387e+12 61768

– Bldg.Type 4 4.4784e+11 4.2773e+12 62012

– Neighborhood 27 1.1457e+12 4.9752e+12 62408

– Overall.Qual 1 1.1960e+12 5.0254e+12 62490

Step: AIC=61692.96

SalePrice ~ Land.Contour + Lot.Config + Land.Slope + Neighborhood +

Condition.1 + Condition.2 + Bldg.Type + House.Style + Overall.Qual +

Overall.Cond + Year.Built + Year.Remod.Add + Roof.Style +

Roof.Matl + Exterior.1st

Df Sum of Sq RSS AIC

3.8313e+12 61693

– Overall.Cond 1 4.5076e+09 3.8358e+12 61694

+ Utilities 2 1.7908e+09 3.8295e+12 61696

– Condition.2 7 3.2008e+10 3.8633e+12 61703

– Land.Contour 3 2.2298e+10 3.8535e+12 61704

– Condition.1 8 4.4153e+10 3.8754e+12 61711

– Lot.Config 4 4.0451e+10 3.8717e+12 61716

– Land.Slope 2 4.3649e+10 3.8749e+12 61722

– Year.Remod.Add 1 6.6441e+10 3.8977e+12 61741

– Year.Built 1 7.0486e+10 3.9017e+12 61744

– House.Style 7 8.8666e+10 3.9199e+12 61746

– Roof.Matl 7 1.0591e+11 3.9372e+12 61759

– Exterior.1st 15 1.3386e+11 3.9651e+12 61764

– Roof.Style 5 1.0995e+11 3.9412e+12 61766

– Bldg.Type 4 4.4838e+11 4.2796e+12 62009

– Neighborhood 27 1.1443e+12 4.9755e+12 62405

– Overall.Qual 1 1.1961e+12 5.0274e+12 62487

# odvajamo ames iz putanje za traženje

detach(ames)

Dijagram raspršenja za varijable.

Dijagram raspršenja varijabli s dodatnim linearnim modelom.

Primijetite da su kategorijske varijable razdijeljene u tzv. dummy varijable, pa se tako svaki nivo faktora pojavljuje kao prediktor u rezultatu.

Koristeći ovih 15 varijabli, objasnili smo 78,84 % varijabiliteta varijable od interesa (cijena nekretnine).

Original Article