Geckaa / APPR-2014-15

Repozitorij z gradivi za predmet Analiza podatkov s programom R v študijskem letu 2014/15
MIT License
0 stars 0 forks source link

Graf #7

Open Geckaa opened 9 years ago

Geckaa commented 9 years ago

Zanima me zakaj mi tak koda za graf: pdf("slike/grafi3.pdf")

prikaz spremembe populacije po letih

leto <- USA$Year pop <- USA$Population.x1000

plot(leto, pop / 1000, xlab = "Year", ylab = "Population x1000000", type = "h", lwd = 3.2, col = "steelblue")

dev.off()

več ne dela, ko pa mi je do zdaj vedno, zdaj pa mi napiše: Error in plot.window(...) : need finite 'xlim' values In addition: Warning messages: 1: In min(x) : no non-missing arguments to min; returning Inf 2: In max(x) : no non-missing arguments to max; returning -Inf 3: In min(x) : no non-missing arguments to min; returning Inf 4: In max(x) : no non-missing arguments to max; returning -In

jaanos commented 9 years ago

Sklepam, da si to razrešil? Meni namreč deluje. Take napake se včasih zgodijo pri risanju v RStudiu, če je okence za graf premajhno.

Geckaa commented 9 years ago

Ja, sem uspel razrešiti. Mislim, da je bilo krivo, da sem počistil delovno okolje in potem pozabil spet zagnati eno zadevo (skratka, butasta napaka na moji strani).

Trenutno me pa zanimata 2 zadevi:

1.) Vsakič ko zaženem svojo vizualizacijo, mi program javi tole:

Checking rgeos availability: FALSE Note: when rgeos is not available, polygon geometry computations in maptools depend on gpclib, which has a restricted licence. It is disabled by default; to enable gpclib, type gpclibPermit() <<

Ko pa izbrišem ta del "gpclibPermit()" in ga ponovno napišem , pa vse dela normalno. SIcer vse deluje potem, je pa to moteče.

2.) Za 4. fazo sem nameraval napovedati funkcijo za rast prebivalstva v ZDA, glede na dane podatke. Program, ki bi narisal graf:

pdf("slike/napoved.pdf") napoved <- function(x,model){predict(model, data.frame(leto=x))} leto <- demographics$Year pop <- demographics$Population.x1000 / 1000 loess <- loess(pop ~ leto) power5 <- lm(pop ~ I(leto^5) + I(leto^4) + I(leto^3) + I(leto^2) + leto) plot(leto, pop, xlim = c(1920, 2030), ylim = c(120, 380), xlab = "Leto", ylab = "Populacija (v mil)") curve(napoved(x, loess), add= TRUE, col = "blue") curve(napoved(x, power5), add = TRUE, col = "red") dev.off() <<

Za moj model (tj. power5) mi javi napako, da bi znalo dati napačne podatke. Točno besedilo ki mi ga javlja, je pa to:

Warning message: In predict.lm(model, data.frame(leto = x)): prediction from a rank-deficient fit may be misleading << Izriše sicer mi, ampak ne vem kaj naj si mislim glede tega.

jaanos commented 9 years ago

Tisto obvestilo se izpiše, ko naložiš knjižnico maptools. Da ga nekoliko zreduciraš, si lahko namestiš še knjižnico rgeos.

Če pogledaš power5$coefficients, lahko vidiš, da sta koeficienta pri linearnem in kvadratnem členu NA. Eden od razlogov, zakaj se to zgodi, je ta, da bi dejanska koeficienta morala biti premajhna. En način, kako to rešiš, je, da premakneš in po potrebi skaliraš os x. V tvojem primeru imaš na osi x leta od 1935 do 2013 - tako lahko odšteješ npr. 1900. Potrebno bo popraviti tudi funkcijo za napoved:

napoved <- function(x,model){predict(model, data.frame(leto=x-1900))}
leto <- demographics$Year - 1900

Seveda bo potrebno pri risanju grafa uporabiti nepremaknjena leta:

plot(demographics$Year, pop, xlim = c(1920, 2030), ylim = c(120, 380), xlab = "Leto", ylab = "Populacija (v mil)")

Krivulje potem rišeš na isti način s popravljeno funkcijo napoved. Sicer pa, glede na to, da je zveza več kot očitno linearna, niti nima smisla poskušati z visokimi potencami (tak popravljen model 5. stopnje bo pokazal, da bi moralo prebivalstvo začeti hitro padati).

Opozoril bi še na komentar, ki sem ga dal pri commitu 6a8bfe039ccc278a71973f966736ff0cee21cb75, zaradi česar se ti program ustavi pri risanju grafov.

Geckaa commented 9 years ago

Spremenil sem imena spremenljivk, nič se več ne konča s piko, upam da zdaj vse dela, na mojem računalniku namreč dela vse. Popravil sem tudi funkcijo napoved.

Rad bi še narisal zemljevid ZDA, ki bi bil pobarvan po tem kako gosto je poseljena posamezna država. V tabelo sem dodal stolpec Pop.Density kjer je gostota preračuana iz podatkov ki sem jih imel na voljo. Zdaj edino ne znam shraniti tega v vektor barv. Probal sem:

gostota <- ZDA[8] barve <- rgb(1, 0, 0, match(states$STATE_NAME, gostota)) Če mi lahko pomagate samo na tem koraku, prosim.

jaanos commented 9 years ago

Če hočeš gostoto predstaviti z enakomerno lestvico (npr. odtenki rdeče), lahko narediš nekaj takega:

gostota <- ZDA[,8]
max.gostota <- max(gostota)
min.gostota <- min(gostota)
norm.gostota <- (gostota - min.gostota)/(max.gostota - min.gostota)
barve <- rgb(1, 0, 0, norm.gostota)
plot(states, col = barve[m])

Pri tem lahko za m uporabiš spremenljivko, ki jo že imaš in si jo dobil z

m <- match(states$STATE_NAME, rownames(ZDA))

Podobno kot pri populaciji tudi tukaj podatki niso enakomerno razporejeni (3. kvartil je šele pri 18,5% največje gostote). Lahko poskusiš nekoliko transformirati lestvico, npr.

barve <- rgb(1, 0, 0, sqrt(norm.gostota))

Seveda dodaj tudi legendo, da bo jasno, kaj katera barva pomeni. Lahko si pa pomagaš tudi z spplot, ki ti sam poračuna barve in nariše lestvico.

Geckaa commented 9 years ago

Legendo sem probal narediti tako:

legend("left", legend = gostota, fill = rgb(1, 0, 0, sqrt(norm.gostota)), cex = 0.5) vendar mi vrne celoten vektor vseh barv, ki si pojavijo na zemljevidu, razporejen po vrstnem redu držav, in polek vsake barve se mi izšpiše vrednost, ki jo ta barva pomeni. To verjamem, da ni vredu. Kako bi naredil eno barvno paleto z samo nekaj vrednostmi?

Imam pa še eno vprašanje. Za 4. fazo sem misliš opraviti še en k-means cluster (recimo da na grafu populacije držav v odvisnosti od gostote prebivalstva). BI bilo to vredu, ali naj kaj drugega?

jaanos commented 9 years ago

Za legendo lahko narediš nekaj takega:

k <- 5 # število stopenj
stopnje <- seq(min.gostota, max.gostota, (max.gostota - min.gostota)/(k-1))
norm.stopnje <- (stopnje - min.gostota)/(max.gostota - min.gostota)
barve.stopnje <- rgb(1, 0, 0, sqrt(norm.stopnje))
legend("left", legend = round(stopnje), fill = barve.stopnje, cex = 0.5)

Seveda se lahko še nekoliko poigraš s samim prikazom legende (namesto round lahko npr. uporabiš signif z ustreznim številom števk) - dodaj tudi pojasnilo, kaj te številke pomenijo (torej enoto - lahko je seveda v naslovu).

k-means clustering bo za 4. fazo v redu - seveda ni smiselno, da ga narediš na enem samem podatku, pač pa lahko uporabiš vse številske podatke, ki jih imaš.

Mimogrede: če hočeš v issueju (ali drugod na GitHubu) prikazati kodo v R-ju, pred njo postavi vrstico

```R

za njo pa



Tako ti bo tudi ustrezno pobarvalo kodo (namesto `R` lahko podaš tudi kak drugi programski jezik, da ti ustrezno obarva, ali pa nič, da ti ne obarva). Če pa hočeš kakšen kos kode vključiti sredi besedila, ga daj med znaka `. S klikom na _Preview_ si lahko ogledaš, kako bo tvoja objava izgledala.
Geckaa commented 9 years ago

Mislim, da sem vglavnem uspel narediti kar je bil moje namen. BI moral še kaj, in ali je to, kar že imam vredu?

jaanos commented 9 years ago

Sedaj vse deluje - ne pozabi pa vključiti 4. faze v glavnem programu.

Še nekaj pripomb imam na izrisane grafe. Pri grafu napovedi populacije bi koristila legenda, iz katere naj bo jasno, katera krivulja je katera - koristno bi bilo tudi izpisati enačbe dobljenih krivulj in pa vsoto kvadratov napake za primerjavo, kateri model se bolje prilega (to je lahko tudi v poročilu). Smiselno bi bilo tudi napovedati, kdaj naj bi po določenem modelu populacija dosegla npr. 400 milijonov (lahko narediš vodoravno črto pri tej vrednosti, pa navpični črti skozi presečišča z vsako krivuljo).

Pri zemljevidih ne pozabi še na naslov (dodaš ga z ukazom title po risanju zemljevida). Smiselno bi bilo narisati še en zemljevid, ki prikazuje delitev na skupine z zadnjih treh grafov. Pri teh ne pozabi še deliti velikosti z npr. 1000 in populacije z 1000000, da se izogneš eksponentni notaciji - pri oseh navedi še uporabljene enote.

Pri razvrščanju v skupine ne pozabi še na normiranje podatkov (funkcija scale) - brez tega bo največjo težo imel podatek z največjimi absolutnimi vrednostmi - pri tebi je to populacija, kot se to lepo vidi z grafov. Svetujem tudi, da ostalih številskih podatkov ne odstraniš, saj bodo v normirani obliki imeli večji vpliv. S pomočjo novih grafov in zemljevida potem to razdelitev tudi interpretiraj - morda se bo pokazala kaka geografska delitev, v vsakem primeru pa lahko delitev interpretiraš na podlagi uporabljenih podatkov.

Geckaa commented 9 years ago

Dopisal sem legende v graf. Prav tako se mi je zdel zanimiv premislek, da bi lahko predvidel kdaj bo populacija dosegla 400 milijonov, tako da sem naredil tudi to. Zemljevidom sem dodal naslov z ukazom title, prav tako dodal en zemljevid, ki države pobarva glede na delitev k-means clusteringa. Imam pa problem z dvema stvarema.

1.) Če bi želel v grafu velikost zdelit z 1000000, mi vrne napako "non-numeric argument to binary operator" Mislil sem, da bi naredil to takole:

plot(data[c("Total.area.in.mi2"/1000, "Population.2013"/1000000)], col = skupina,
     xlab = "Velikost", ylab = "Populacija",
     main = "Grupiranje - Odvisnot populacije od velikosti države")

2.) Če zemljevid narišem in pobarvam z ukazom:

plot(states, col = skupina)

sam ne morem vplivat na barve, s katerimi pobarva, in ker imam delitev na 3 skupine mi vedno večino pobarva s črno barvo, kar stvar naredi nepregledno. Kako to rešiti?

jaanos commented 9 years ago

1.) Tukaj poskušaš deliti (imenske) indekse, kar seveda ne gre. Treba bo deliti podatke:

plot(data[,"Total.area.in.mi2"]/1000, data[,"Population.2013"]/1000000, col = skupina,
     xlab = "Velikost", ylab = "Populacija",
     main = "Grupiranje - Odvisnot populacije od velikosti države")

2.) Treba bo definirati vektor barv in ga indeksirati z vektorjem skupina. Takole bi šlo:

barve <- c("red", "green", "blue", "yellow")
plot(states, col = barve[skupina])

Na enak način potem določiš barve tudi pri grafih.

Geckaa commented 9 years ago

Najlepša hvala, to sem zdaj popravil. Moram dodati oziroma popraviti še kaj, ali je vredu to kar imam?

jaanos commented 9 years ago

Izgleda v redu - v poročilo bi dodal edino še enačbi krivulj za vsak model.

Če želiš pristopiti k zagovoru, odpri issue na repozitoriju za zagovore s šifro commita, s profesorjem se pa domeni za termin zagovora.

Geckaa commented 9 years ago

Odprl sem issue, lahko samo preverite če je vredu?

jaanos commented 9 years ago

Je v redu, sem dodal.