+ Näytä koodi- Piilota koodi
library(sorvi)
library(OpasnetBaseUtils)
library(OpasnetUtils)
library(xtable)
library(ggplot2)
kunnat <- data.frame(Kunta = c("Akaa", "Alajärvi", "Alavieska", "Alavus", "Asikkala", "Askola", "Aura", "Brändö", "Eckerö",
"Enonkoski", "Enontekiö", "Espoo", "Eura", "Eurajoki", "Evijärvi", "Finström", "Forssa", "Föglö",
"Geta", "Haapajärvi", "Haapavesi", "Hailuoto", "Halsua", "Hamina", "Hammarland", "Hankasalmi",
"Hanko", "Harjavalta", "Hartola", "Hattula", "Haukipudas", "Hausjärvi", "Heinävesi", "Helsinki",
"Vantaa", "Hirvensalmi", "Hollola", "Honkajoki", "Huittinen", "Humppila", "Hyrynsalmi",
"Hyvinkää", "Hämeenkyrö", "Hämeenlinna", "Heinola", "Ii", "Iisalmi", "Iitti", "Ikaalinen",
"Ilmajoki", "Ilomantsi", "Inari", "Inkoo", "Isojoki", "Isokyrö", "Imatra", "Jalasjärvi",
"Janakkala", "Joensuu", "Jokioinen", "Jomala", "Joroinen", "Joutsa", "Juankoski", "Juuka",
"Juupajoki", "Juva", "Jyväskylä", "Jämijärvi", "Jämsä", "Järvenpää", "Kaarina", "Kaavi",
"Kajaani", "Kalajoki", "Kangasala", "Kangasniemi", "Kankaanpää", "Kannonkoski", "Kannus",
"Karijoki", "Karjalohja", "Karkkila", "Karstula", "Karvia", "Kaskinen", "Kauhajoki", "Kauhava",
"Kauniainen", "Kaustinen", "Keitele", "Kemi", "Keminmaa", "Kemiönsaari", "Kempele", "Kerava",
"Kerimäki", "Kesälahti", "Keuruu", "Kihniö", "Kiikoinen", "Kiiminki", "Kinnula", "Kirkkonummi",
"Kitee", "Kittilä", "Kiuruvesi", "Kivijärvi", "Kokemäki", "Kokkola", "Kolari", "Konnevesi",
"Kontiolahti", "Korsnäs", "Hämeenkoski", "Koski Tl", "Kotka", "Kouvola", "Kristiinankaupunki",
"Kruunupyy", "Kuhmo", "Kuhmoinen", "Kumlinge", "Kuopio", "Kuortane", "Kurikka", "Kustavi",
"Kuusamo", "Outokumpu", "Kyyjärvi", "Kärkölä", "Kärsämäki", "Kökar", "Köyliö", "Kemijärvi",
"Lahti", "Laihia", "Laitila", "Lapinlahti", "Lappajärvi", "Lappeenranta", "Lapinjärvi", "Lapua",
"Laukaa", "Lavia", "Lemi", "Lemland", "Lempäälä", "Leppävirta", "Lestijärvi", "Lieksa", "Lieto",
"Liminka", "Liperi", "Loimaa", "Loppi", "Loviisa", "Luhanka", "Lumijoki", "Lumparland", "Larsmo",
"Luumäki", "Luvia", "Lohja", "Länsi-Turunmaa", "Maalahti", "Maaninka", "Mariehamn", "Marttila",
"Masku", "Merijärvi", "Merikarvia", "Miehikkälä", "Mikkeli", "Muhos", "Multia", "Muonio",
"Mustasaari", "Muurame", "Mynämäki", "Myrskylä", "Mäntsälä", "Mänttä-Vilppula", "Mäntyharju",
"Naantali", "Nakkila", "Nastola", "Nilsiä", "Nivala", "Nokia", "Nousiainen", "Nummi-Pusula",
"Nurmes", "Nurmijärvi", "Närpiö", "Orimattila", "Oripää", "Orivesi", "Oulainen", "Oulu",
"Oulunsalo", "Padasjoki", "Paimio", "Paltamo", "Parikkala", "Parkano", "Pelkosenniemi", "Perho",
"Pertunmaa", "Petäjävesi", "Pieksämäki", "Pielavesi", "Pietarsaari", "Pedersören kunta",
"Pihtipudas", "Pirkkala", "Polvijärvi", "Pomarkku", "Pori", "Pornainen", "Posio", "Pudasjärvi",
"Pukkila", "Punkaharju", "Punkalaidun", "Puolanka", "Puumala", "Pyhtää", "Pyhäjoki", "Pyhäjärvi",
"Pyhäntä", "Pyhäranta", "Pälkäne", "Pöytyä", "Porvoo", "Raahe", "Raasepori", "Raisio",
"Rantasalmi", "Ranua", "Rauma", "Rautalampi", "Rautavaara", "Rautjärvi", "Reisjärvi", "Riihimäki",
"Ristiina", "Ristijärvi", "Rovaniemi", "Ruokolahti", "Ruovesi", "Rusko", "Rääkkylä", "Saarijärvi",
"Salla", "Salo", "Saltvik", "Sastamala", "Sauvo", "Savitaipale", "Savonlinna", "Savukoski",
"Seinäjoki", "Sievi", "Siikainen", "Siikajoki", "Siikalatva", "Siilinjärvi", "Simo", "Sipoo",
"Siuntio", "Sodankylä", "Soini", "Somero", "Sonkajärvi", "Sotkamo", "Sottunga", "Sulkava", "Sund",
"Suomenniemi", "Suomussalmi", "Suonenjoki", "Sysmä", "Säkylä", "Vaala", "Taipalsaari",
"Taivalkoski", "Taivassalo", "Tammela", "Tampere", "Tarvasjoki", "Tervo", "Tervola", "Teuva",
"Tohmajärvi", "Toholampi", "Toivakka", "Tornio", "Turku", "Pello", "Tuusniemi", "Tuusula",
"Tyrnävä", "Töysä", "Ulvila", "Urjala", "Utajärvi", "Utsjoki", "Uurainen", "Uusikaarlepyy",
"Uusikaupunki", "Vaasa", "Valkeakoski", "Valtimo", "Varkaus", "Vehmaa", "Vesanto", "Vesilahti",
"Veteli", "Vieremä", "Vihanti", "Vihti", "Viitasaari", "Vimpeli", "Virolahti", "Virrat", "Vårdö",
"Vähäkyrö", "Vöyri", "Yli-Ii", "Ylitornio", "Ylivieska", "Ylöjärvi", "Ypäjä", "Ähtäri",
"Äänekoski"),
Maakunta = c("Pirkanmaa", "Etelä-Pohjanmaa", "Pohjois-Pohjanmaa", "Etelä-Pohjanmaa",
"Päijät-Häme", "Itä-Uusimaa", "Varsinais-Suomi", "Ahvenanmaa", "Ahvenanmaa", "Etelä-Savo",
"Lappi", "Uusimaa", "Satakunta", "Satakunta", "Etelä-Pohjanmaa", "Ahvenanmaa", "Kanta-Häme",
"Ahvenanmaa", "Ahvenanmaa", "Pohjois-Pohjanmaa", "Pohjois-Pohjanmaa", "Pohjois-Pohjanmaa",
"Keski-Pohjanmaa", "Kymenlaakso", "Ahvenanmaa", "Keski-Suomi", "Uusimaa", "Satakunta",
"Päijät-Häme", "Kanta-Häme", "Pohjois-Pohjanmaa", "Kanta-Häme", "Etelä-Savo", "Uusimaa",
"Uusimaa", "Etelä-Savo", "Päijät-Häme", "Satakunta", "Satakunta", "Kanta-Häme", "Kainuu",
"Uusimaa", "Pirkanmaa", "Kanta-Häme", "Päijät-Häme", "Pohjois-Pohjanmaa", "Pohjois-Savo",
"Kymenlaakso", "Pirkanmaa", "Etelä-Pohjanmaa", "Pohjois-Karjala", "Lappi", "Uusimaa",
"Etelä-Pohjanmaa", "Pohjanmaa", "Etelä-Karjala", "Etelä-Pohjanmaa", "Kanta-Häme",
"Pohjois-Karjala", "Kanta-Häme", "Ahvenanmaa", "Etelä-Savo", "Keski-Suomi", "Pohjois-Savo",
"Pohjois-Karjala", "Pirkanmaa", "Etelä-Savo", "Keski-Suomi", "Satakunta", "Keski-Suomi",
"Uusimaa", "Varsinais-Suomi", "Pohjois-Savo", "Kainuu", "Pohjois-Pohjanmaa", "Pirkanmaa",
"Etelä-Savo", "Satakunta", "Keski-Suomi", "Keski-Pohjanmaa", "Etelä-Pohjanmaa", "Uusimaa",
"Uusimaa", "Keski-Suomi", "Satakunta", "Pohjanmaa", "Etelä-Pohjanmaa", "Etelä-Pohjanmaa",
"Uusimaa", "Keski-Pohjanmaa", "Pohjois-Savo", "Lappi", "Lappi", "Varsinais-Suomi",
"Pohjois-Pohjanmaa", "Uusimaa", "Etelä-Savo", "Pohjois-Karjala", "Keski-Suomi", "Pirkanmaa",
"Satakunta", "Pohjois-Pohjanmaa", "Keski-Pohjanmaa", "Uusimaa", "Pohjois-Karjala", "Lappi",
"Pohjois-Savo", "Keski-Suomi", "Satakunta", "Keski-Pohjanmaa", "Lappi", "Keski-Suomi",
"Pohjois-Karjala", "Pohjanmaa", "Päijät-Häme", "Varsinais-Suomi", "Kymenlaakso", "Kymenlaakso",
"Pohjanmaa", "Pohjanmaa", "Kainuu", "Keski-Suomi", "Ahvenanmaa", "Pohjois-Savo",
"Etelä-Pohjanmaa", "Etelä-Pohjanmaa", "Varsinais-Suomi", "Pohjois-Pohjanmaa", "Pohjois-Karjala",
"Keski-Suomi", "Päijät-Häme", "Pohjois-Pohjanmaa", "Ahvenanmaa", "Satakunta", "Lappi",
"Päijät-Häme", "Pohjanmaa", "Varsinais-Suomi", "Pohjois-Savo", "Etelä-Pohjanmaa", "Etelä-Karjala",
"Itä-Uusimaa", "Etelä-Pohjanmaa", "Keski-Suomi", "Satakunta", "Etelä-Karjala", "Ahvenanmaa",
"Pirkanmaa", "Pohjois-Savo", "Keski-Pohjanmaa", "Pohjois-Karjala", "Varsinais-Suomi",
"Pohjois-Pohjanmaa", "Pohjois-Karjala", "Varsinais-Suomi", "Kanta-Häme", "Itä-Uusimaa",
"Keski-Suomi", "Pohjois-Pohjanmaa", "Ahvenanmaa", "Pohjanmaa", "Etelä-Karjala", "Satakunta",
"Uusimaa", "Varsinais-Suomi", "Pohjanmaa", "Pohjois-Savo", "Ahvenanmaa", "Varsinais-Suomi",
"Varsinais-Suomi", "Pohjois-Pohjanmaa", "Satakunta", "Kymenlaakso", "Etelä-Savo",
"Pohjois-Pohjanmaa", "Keski-Suomi", "Lappi", "Pohjanmaa", "Keski-Suomi", "Varsinais-Suomi",
"Itä-Uusimaa", "Itä-Uusimaa", "Pirkanmaa", "Etelä-Savo", "Varsinais-Suomi", "Satakunta",
"Päijät-Häme", "Pohjois-Savo", "Pohjois-Pohjanmaa", "Pirkanmaa", "Varsinais-Suomi", "Uusimaa",
"Pohjois-Karjala", "Uusimaa", "Pohjanmaa", "Päijät-Häme", "Varsinais-Suomi", "Pirkanmaa",
"Pohjois-Pohjanmaa", "Pohjois-Pohjanmaa", "Pohjois-Pohjanmaa", "Päijät-Häme", "Varsinais-Suomi",
"Kainuu", "Etelä-Karjala", "Pirkanmaa", "Lappi", "Keski-Pohjanmaa", "Etelä-Savo", "Keski-Suomi",
"Etelä-Savo", "Pohjois-Savo", "Pohjanmaa", "Pohjanmaa", "Keski-Suomi", "Pirkanmaa",
"Pohjois-Karjala", "Satakunta", "Satakunta", "Itä-Uusimaa", "Lappi", "Pohjois-Pohjanmaa",
"Itä-Uusimaa", "Etelä-Savo", "Pirkanmaa", "Kainuu", "Etelä-Savo", "Kymenlaakso",
"Pohjois-Pohjanmaa", "Pohjois-Pohjanmaa", "Pohjois-Pohjanmaa", "Varsinais-Suomi", "Pirkanmaa",
"Varsinais-Suomi", "Itä-Uusimaa", "Pohjois-Pohjanmaa", "Uusimaa", "Varsinais-Suomi", "Etelä-Savo",
"Lappi", "Satakunta", "Pohjois-Savo", "Pohjois-Savo", "Etelä-Karjala",
"Keski-Pohjanmaa/Pohjois-Pohjanmaa", "Kanta-Häme", "Etelä-Savo", "Kainuu", "Lappi",
"Etelä-Karjala", "Pirkanmaa", "Varsinais-Suomi", "Pohjois-Karjala", "Keski-Suomi", "Lappi",
"Varsinais-Suomi", "Ahvenanmaa", "Pirkanmaa", "Varsinais-Suomi", "Etelä-Karjala", "Etelä-Savo",
"Lappi", "Etelä-Pohjanmaa", "Keski-Pohjanmaa", "Satakunta", "Pohjois-Pohjanmaa",
"Pohjois-Pohjanmaa", "Pohjois-Savo", "Lappi", "Itä-Uusimaa", "Uusimaa", "Lappi",
"Etelä-Pohjanmaa", "Varsinais-Suomi", "Pohjois-Savo", "Kainuu", "Ahvenanmaa", "Etelä-Savo",
"Ahvenanmaa", "Etelä-Karjala", "Kainuu", "Pohjois-Savo", "Päijät-Häme", "Satakunta", "Kainuu",
"Etelä-Karjala", "Pohjois-Pohjanmaa", "Varsinais-Suomi", "Kanta-Häme", "Pirkanmaa",
"Varsinais-Suomi", "Pohjois-Savo", "Lappi", "Etelä-Pohjanmaa", "Pohjois-Karjala",
"Keski-Pohjanmaa", "Keski-Suomi", "Lappi", "Varsinais-Suomi", "Lappi", "Pohjois-Savo", "Uusimaa",
"Pohjois-Pohjanmaa", "Etelä-Pohjanmaa", "Satakunta", "Pirkanmaa", "Pohjois-Pohjanmaa", "Lappi",
"Keski-Suomi", "Pohjanmaa", "Varsinais-Suomi", "Pohjanmaa", "Pirkanmaa", "Pohjois-Karjala",
"Pohjois-Savo", "Varsinais-Suomi", "Pohjois-Savo", "Pirkanmaa", "Keski-Pohjanmaa", "Pohjois-Savo",
"Pohjois-Pohjanmaa", "Uusimaa", "Keski-Suomi", "Etelä-Pohjanmaa", "Kymenlaakso", "Pirkanmaa",
"Ahvenanmaa", "Pohjanmaa", "Pohjanmaa", "Pohjois-Pohjanmaa", "Lappi", "Pohjois-Pohjanmaa",
"Pirkanmaa", "Kanta-Häme", "Etelä-Pohjanmaa", "Keski-Suomi"))
# print(xtable(kunnat), type = 'html')
#head(kunnat)
# Luonnos uudemmaksi tavaksi kaivaa data.
rad <- op_baseGetData("opasnet_base", "Op_fi2759")
colnames(rad)[colnames(rad) == "Havainto"] <- "Observation"
rad <- tidy(rad, objname = "radon", direction = "long")
print(xtable(rad), type = 'html')
radon <- op_baseGetData("opasnet_base", "op_fi2759")[, 3:6]
radon$Result <- ifelse(is.na(radon$Result.Text), radon$Result, as.numeric(as.character(radon$Result.Text)))
radon <- radon[, 1:3]
colnames(radon)[3] <- "r"
#head(radon)
radon <- radon[substr(radon$Havainto, 1, 7) != "kaikki.", ]
radon <- radon[radon$Maakunta != "Yhteensä", ]
#print(xtable(reshape(radon, idvar = "Maakunta", timevar = "Havainto", direction = "wide")), type = 'html')
radon$Talo <- ifelse(substr(radon$Havainto, 1, 3) == "pt.", "Pientalo", "Kerrostalo")
radon$P <- substr(radon$Havainto, 4, 10)
radon <- radon[, colnames(radon)!="Havainto"]
radon <- reshape(radon, idvar = c("Maakunta", "Talo"), timevar = "P", direction = "wide")
radon[is.na(radon)] <- 0
#print(xtable(radon), type = 'html')
j <- 2/3 # c( rep(2/3, 8), 3/4, rep(2/3, 9), 3/4, rep(2/3, 21))
X2 <- radon$r.max
for(i in 1:nrow(radon)) {X2[i] <- min(radon$r.max[i], (radon$r.med[i] + 201)/2)}
X <- data.frame(P1 = radon$r.med, P2 = X2, P3 = radon$r.max, P4 = 201, P5 = 351, P6 = 401, P7 = 601, P8 = 801)
P <- data.frame(
P1 = 50,
P2 = (50 - radon$r.y200),
P3 = (radon$r.y800)*(1-j),
P4 = (radon$r.y200 - radon$r.y400)*j,
P5 = (radon$r.y200 - radon$r.y400)*(1-j),
P6 = (radon$r.y400 - radon$r.y800)*j,
P7 = (radon$r.y400 - radon$r.y800)*(1-j),
P8 = (radon$r.y800)*j)/100
#X <- data.frame(P1 = radon$r.Mediaani, P2 = radon$r.Maksimi, P3 = 201, P4 = 351, P5 = 401, P6 = 601, P7 = 801)
#P <- data.frame(P1 = 50, P2 = radon$r.Yli800*(1-j), P3 = radon$r.Yli200*j, P4 = radon$r.Yli200*(1-j), P5 = radon$r.Yli400*j, P6 = radon$r.Yli400*(1-j), P7 = radon$r.Yli800*j)/100
Xsort <- X
for(i in 1:nrow(X)){Xsort[i, ] <- sort.int(t(X[i, ]), index.return=TRUE)[[2]]}
for(i in 1:nrow(X)){X[i, ] <- X[i, t(Xsort[i, ])]}
for(i in 1:nrow(X)){P[i, ] <- P[i, t(Xsort[i, ])]}
#head(radon)
#print(xtable(Xsort), type = 'html')
#print(xtable(P), type = 'html')
#head(P)
X$P1 <- (radon$r.ka - rowSums((X*P)[, -1])) / 0.5
#print(xtable(X), type = 'html')
#rowSums(X*P)
#rowSums(P)
intP <- P
for(i in 2:ncol(P)){intP[, i] <- rowSums(P[, 1:i])}
#intP
out <- data.frame()
for(j in 1:n){
radoniter <- intP[, 1]
iter <- rowSums(runif(nrow(intP), 0, 1) > intP)+1
for(i in 1:nrow(intP)){radoniter[i] <- X[i, iter[i]]}
radon2 <- data.frame(Maakunta = radon$Maakunta, Talo = radon$Talo, obs = n, Radon = radoniter)
radon2
out <- rbind(out, radon2)
}
#radon2
#head(out)
talot <- op_baseGetData("opasnet_base", "Op_fi2761")[, -c(1, 2, 8)]
talot <- talot[talot$Talotyyppi != "Kaikki talotyypit" & talot$Kunta !="KOKO MAA - HELA LANDET" & talot$Asuntokunta != "Kaikki asuntokunnat", ]
talot$Talotyyppi <- ifelse(talot$Talotyyppi=="Asuinkerrostalo", "Kerrostalo", "Pientalo")
#head(talot)
väki <- data.frame(Asuntokunta = c("1 henk.", "2 henk.", "3 henk.", "4 henk.", "5 henk.", "6 henk.", "7+ henk."), n = 1:7)
talot <- merge(talot, väki)
#head(talot)
talot$Result <- talot$Result * talot$n
talot <- as.data.frame(as.table(tapply(talot$Result, list(talot$Kunta, talot$Talotyyppi), sum)))
colnames(talot) <- c("Kunta", "Talo", "Asukkaita")
#head(talot)
talot <- merge(talot, kunnat)
#head(talot)
out <- merge(out, talot)
#head(out)
rr.radon <- 1.0016
rr.tupakka <- 20
tausta.altistus <- 5
tausta.sairastuvuus <- data.frame( Maakunta = rep(c("Uusimaa", "Itä-Uusimaa", "Varsinais-Suomi", "Satakunta", "Kanta-Häme",
"Pirkanmaa", "Päijät-Häme", "Kymenlaakso", "Etelä-Karjala", "Etelä-Savo", "Itä-Savo", "Pohjois-Karjala", "Pohjois-Savo",
"Keski-Suomi", "Pohjanmaa", "Vaasa", "Keski-Pohjanmaa", "Pohjois-Pohjanmaa", "Kainuu", "Länsi-Pohja", "Lappi", "Ahvenanmaa", "Koko maa"), 2),
Sukupuoli = (rep( c("Miehet", "Naiset"), each = 23)), Sairastuvuus = c(29.5, 29.5, 26.3, 30.2, 30.3, 26.5, 26.9, 26.6, 29.0, 26.1, 28.7, 29.0, 26.7, 24.7, 26.1, 26.1, 29.8, 30.9, 28.3, 33.7, 34.3, 27.9, 28.2, 11.8, 11.8, 9.2, 7.7, 9.0, 7.7, 7.7, 9.5, 7.2, 7.4, 6.6, 5.8, 5.4, 6.0, 5.9, 8.7, 9.3, 8.7, 7.9, 10.1, 11.3, 13.9, 8.8))
tupakoivat <- data.frame(Sukupuoli = rep(c("Miehet", "Naiset"), 2), Tupakoivien.osuus = rep(c(0.22, 0.16), 2), Tupakka = rep(c("Kyllä", "Ei"), each = 2))
out <- merge(out, tausta.sairastuvuus)
out <- merge(out, tupakoivat)
out$Asukkaita <- out$Asukkaita/2 * ifelse(out$Tupakka == "Kyllä", out$Tupakoivien.osuus, 1-out$Tupakoivien.osuus)
out$Sairastuvuus <- out$Sairastuvuus * (1 - out$Tupakoivien.osuus * (rr.tupakka - 1) / (out$Tupakoivien.osuus * (rr.tupakka - 1) +1))
out$Kokonaisvaikutus <- exp(log(rr.tupakka)*ifelse(out$Tupakka == "Kyllä", 1,0) + log(rr.radon)*out$Radon)
out$Radonvaikutus <- out$Kokonaisvaikutus - exp(log(rr.tupakka)*ifelse(out$Tupakka == "Kyllä", 1,0))
out$Radonvaikutus <- out$Radonvaikutus * out$Sairastuvuus/100000 * out$Asukkaita
out$Kokonaisvaikutus <- out$Kokonaisvaikutus * out$Sairastuvuus/100000 * out$Asukkaita
print(xtable(out[1:20, ]), type = 'html', html.table.attributes="border=1 class='sortable'")
print("Yksikköriski: tapauksia tuhatta asukasta kohti vuodessa")
temp <- as.data.frame(as.table(tapply(out$Radonvaikutus/out$Asukkaita, list(out$Maakunta, out$Tupakka, out$Talo), sum)/n))
temp <- temp[!is.na(temp$Freq), ]
#head(temp)
temp$Freq <- temp$Freq*1000
colnames(temp) <- c("Maakunta", "Tupakka", "Talotyyppi", "Riski")
print(xtable(temp), type= 'html', html.table.attributes="border=1 class='sortable'")
print("Radonin aiheuttama lisäriski, tapausta vuodessa")
print(xtable(as.data.frame(as.table(((tapply(out$Radonvaikutus, list(out$Talo, out$Tupakka), sum)/n))))), type = 'html', html.table.attributes="border=1 class='sortable'")
ggplot(out, aes(x = Maakunta, weight = Kokonaisvaikutus, fill = Tupakka)) + geom_bar(position = "stack") # position = "dodge"
out$Radonvaikutus <- ifelse(out$Radonvaikutus < 0, 0, out$Radonvaikutus)
out$Radonsata <- round(as.numeric(as.character(out$Radon)), -2)
ggplot(out, aes(x = Radon, weight = Radonvaikutus)) + geom_density(alpha = 1)
ggplot(out, aes(x = Radonsata, weight = Radonvaikutus, fill = Talo)) + geom_bar(position = "stack")
###########################################################
# (C) 2011 Leo Lahti <leo.lahti[at]iki.fi> All rights reserved.
# License: FreeBSD, http://en.wikipedia.org/wiki/BSD_licenses
# Tama esimerkki on testattu sorvi-paketin versiolla 0.1.23
# Esimerkki Suomen kuntatason vaestonkasvutilastojen (Tilastokeskus)
# visualisoinnista Maanmittauslaitoksen karttadatalla (vuonna 2010)
###############################################
# Lue Suomen kuntarajat SpatialPolygon-muodossa
# (C) Maanmittauslaitos 2011
# http://www.maanmittauslaitos.fi/aineistot-palvelut/digitaaliset-tuotteet/ilmaiset-aineistot/hankinta
#LoadData(MML)
#sp <- MML[["1_milj_Shape_etrs_shape"]][["kunta1_p"]]
sp <- LoadMML(data.id = "kunta1_p", resolution = "1_milj_Shape_etrs_shape")
#################################################
# Lue kuntatason vaestonkasvutiedot tilastokeskuksen StatFin-tietokannasta
# http://www.stat.fi/tup/statfin/index.html
# PC Axis-muodossa
px <- statfin.px("vrm/synt/080_synt_tau_203.px")
# statfin.px-funktio käyttää soRvin read.px-funktiota hakemaan dataa URLista, joka alkaa "http://pxweb2.stat.fi/database/StatFin/"
# URLin loppuosa pitää antaa statfin.px:lle parametrina.
# Poimi taulukosta halutut tiedot
# pxs <- subset(as.data.frame(px), Väestönmuutos.ja.väkiluku == "Luonnollinen väestönlisäys" & Vuosi == year)
# pxs <- subset(as.data.frame(px), Väestönmuutos.ja.väkiluku == parameter & Vuosi == year)
# Putsaa data
# vaestonkasvu <- preprocess.px(pxs)
# head(vaestonkasvu)
################################################
vaestonkasvu <- out[out$Talo == talo & out$Tupakka == tupakka, ]
vaestonkasvu <- as.data.frame(as.table(tapply(vaestonkasvu[, vaste], list(vaestonkasvu$Kunta), sum)/n))
colnames(vaestonkasvu) <- c("Kunta", "Vaste")
head(vaestonkasvu)
vaestonkasvu[is.na(vaestonkasvu)] <- 0
##############################################
# Lisaa tiedot karttaobjektiin
sp@data$vaestonkasvu <- vaestonkasvu$Vaste[match(sp$Kunta.FI, vaestonkasvu$Kunta)]
# Korvaa puuttuvat arvot nollalla
sp[["vaestonkasvu"]][is.na(sp[["vaestonkasvu"]])] <- 0
#colnames(out)[vaste]
#vaestonkasvu$Vaste
################################################
# Maarittele varipaletti
my.palette <- colorRampPalette(c("white", "red"), space = "rgb") # , "blue"
ncol <- 40 # Number of colors
#################################################
# Piirra kuva
varname <- "vaestonkasvu"
int <- max(abs(vaestonkasvu$Vaste))
int
q <- spplot(sp, varname,
col.regions = my.palette(ncol),
main = paste(colnames(out)[vaste], talo, "Tupakointi: ", tupakka, sep = ", "),
colorkey = TRUE,
lwd = .4,
col = "black",
at = seq(0, 0 + int, length = ncol)
)
print(q)
| |