+ Näytä koodi- Piilota koodi
library(OpasnetUtils)
library(ggplot2)
objects.latest('Op_en6007', code_name = 'answer') # [[OpasnetUtils/Drafts]] fetches fillna and collapsemarg function.
###### OSA 1: VÄESTÖ
### TALOT
objects.latest('Op_fi2761', code_name = 'alusta') # [[Talotyypit Suomessa]], ovariable talot.
# Yksinkertaista talotyyppiluokittelua ja summaa taloista pois perhekoko.
#talot@output <- talot@output[talot@output$Kunta == "Akaa" , ] ########### POISTA TÄMÄ
levels(talot@output$Talotyyppi) <- ifelse(
levels(talot@output$Talotyyppi) == "Asuinkerrostalo",
"Kerrostalo",
"Pientalo"
)
talot <- collapsemarg(talot, "Vuosi", fun = "pick", picks = list(c("2010")))
talot <- collapsemarg(talot, c("Asuntokunta"), "sum")
### KUNNAT
objects.latest('Op_fi3907', code_name = 'alusta') # [[Suomen kunnat]], ovariable kunnat.
kunnat <- EvalOutput(kunnat, N = N)
### TUPAKOINTI JA ALUEVASTAAVUUS
tupakointi <- EvalOutput(Ovariable("tupakointi", # Tupakoivien osuus väestöstä. Result-sarake summautuu ykköseen.
data = data.frame(
Sukupuoli = rep(c("Miehet", "Naiset"), 2),
Tupakka = rep(c("Kyllä", "Ei"), each = 2),
Result = c(0.11, 0.08, 0.39, 0.42)
)
), N = N)
aluevastaavuus <- EvalOutput(Ovariable("aluevastaavuus",
data = {
dat <- tidy(opbase.data("Op_fi2760.maakunnat_ja_sairaanhoitopiirit"))
colnames(dat)[colnames(dat) == "Result"] <- "Sairaanhoitopiiri"
dat$Result <- 1
dat
}
))
### VÄESTÖN LASKENTA
väestö <- kunnat * talot * aluevastaavuus * tupakointi
if(verbose) {
cat("Väestö\n")
oprint(head(väestö@output))
cat("Kunnat\n")
oprint(head(kunnat@output))
cat("Talot\n")
oprint(head(talot@output))
cat("Aluevastaavuus\n")
oprint(head(aluevastaavuus@output))
cat("Tupakointi\n")
oprint(head(tupakointi@output))
}
###### OSA 2: RISKISUHTEET
### ALTISTUS
# RADON
objects.latest('Op_fi2759', code_name = 'alusta') # [[Radon sisäilmassa]], ovariable radonpit, data.frame radon.
radonpit <- EvalOutput(radonpit, N = N)
#radonpit@output <- radonpit@output[radonpit@output$Maakunta == "Pirkanmaa" , ] ############# POISTA TÄMÄ
# TUPAKKA. Väestö on jo jaettu tupakoiviin ja tupakoimattomiin.
tupakkaaltistus <- EvalOutput(Ovariable("tupakkaaltistus",
data = data.frame(Tupakka = c("Kyllä", "Ei"), Result = c(1, 0))
), N = N)
# Altistuksessa lisätään uusi sarake Altiste kuvaamaan sitä, mikä altisteen altistumisesta tai riskistä on kyse ko. rivillä.
altistus <- orbind(
tupakkaaltistus * EvalOutput(Ovariable(data = data.frame(Altiste = "Tupakka", Result = 1)), N = N),
radonpit * EvalOutput(Ovariable(data = data.frame(Altiste = "Radon", Result = 1)), N = N)
)
altistus <- EvalOutput(Ovariable("altistus", data = altistus), N = N)
altistus@output <- fillna(altistus@output, c("Tupakka", "Maakunta", "Talotyyppi"))
if(verbose) {
cat("Altistus\n")
oprint(head(altistus@output))
}
### SUHTEELLINEN RISKI
# Annosvasteiden perustiedot.
ERF <- EvalOutput(Ovariable("ERF", data = data.frame(Altiste = c("Radon", "Tupakka"), Result = c(1.0016, 20))), N = N)
tausta.altistus <- 0 # Radonille voisi olla 5.
RR <- exp(log(ERF) * (altistus - tausta.altistus)) # Relative risk given the exposures
if(verbose) {
cat("ERF\n")
oprint(head(ERF@output))
cat("RR\n")
oprint(head(RR@output))
}
### SYÖPÄKUOLLEISUUS
objects.latest('Op_fi3912', code_name = 'alusta') # [[Syöpäkuolleisuus Suomessa]], ovariable syopakuolleisuus.
syopakuolleisuus@data <- syopakuolleisuus@data[syopakuolleisuus@data$Primaaripaikka == "Keuhkot, henkitorvi" , ] # ICD.10.koodi == "C33-34"
syopakuolleisuus <- EvalOutput(syopakuolleisuus, N = N)
if(verbose) {
cat("Syöpäkuolleisuus\n")
oprint(head(syopakuolleisuus@output))
}
### TAUSTARISKI
# Calculate subgroup-specific background risks by first calculating a temporary burden estimate.
# temp is a population-weighted average of relative risks. The actual disease risk is divided by this value and then multiplied
# by the subgroup-specific relative risk. In this way, the total burden of all subgroups equals to what is actually seen.
temp <- collapsemarg(RR, c("Altiste"), prod)
temp1 <- temp * väestö # Population-weighted sum of the relative risk.
#KORJAA temp1 <- collapsemarg(temp1, c("Maakunta"Sairaanhoitopiiri", "Sukupuoli", "Iter"), sum) # Aggregate to the same indices as kokonaisvaikutus.
temp2 <- (temp * 0 + 1) * väestö # Population-weighted sum of ones with the same dimensions.
#KORJAA temp2 <- collapsemarg(temp2, c("Sairaanhoitopiiri", "Sukupuoli", "Iter"), sum)
tausta <- syopakuolleisuus / (temp1 / temp2) # syopakuolleisuus without any exposure
if(verbose) {
cat("Taustariski\n")
oprint(head(tausta@output))
}
###### OSA 3: TAUTIKUORMA
### ARVIO PERUSTUEN POPULATION ATTRIBUTABLE FRACTIONIIN
# Tässä oletetaan, että riskisuhteiden avulla voidaan suoraan laskea eri teijöiden osuus tautikuormasta.
# Luku on yliarvio, koska jokainen altiste lasketaan ikään kuin muita altisteita ei olisi tautia aiheuttamassa.
# Etuna on, että jokaiselle altisteelle saadaan arvio erikseen.
väestöosuus <- 1
PAF <- väestöosuus * (RR - 1) / (väestöosuus * (RR - 1) + 1)
vaikutuspaf <- väestö * syopakuolleisuus * PAF
# Onko tarpeen integroida yli muiden indeksien?
#vaikutuspaf <- collapsemarg(vaikutuspaf, c("Iter", "Kunta, ", "Tupakka", "Talotyyppi", "Altiste"), sum)
### ARVIO PERUSTUEN TAUSTARISKIIN
# Tässä oletetaan, että on olemassa populaatiokohtainen taustariski, jota ei altistumisen takia havaita.
# Se voidaan kuitenkin laskea ottamalla havaittu tautikuorma ja poistamalla altistumisen vaikutus.
# Arvio lienee tarkempi kuin PAF-menetelmällä, mutta altistekohtainen tieto menetetään.
RRkok <- collapsemarg(RR, c("Iter", "Tupakka", "Maakunta", "Talotyyppi"), prod)
vaikutusrr <- väestö * tausta * RRkok
#vaikutusrr <- collapsemarg(vaikutusrr, c("Iter", "Kunta", "Tupakka", "Talotyyppi"), sum)
if(verbose) {
cat("PAF\n")
oprint(head(PAF@output))
cat("VaikutusPAF\n")
oprint(head(vaikutus@output))
cat("Kokonaisvaikutus PAF-menetelmällä\n")
oprint(head(kokonaisvaikutuspaf@output))
cat("RRkok\n")
oprint(head(RRkok@output))
cat("VaikutusRR\n")
oprint(head(vaikutusrr@output))
cat("Kokonaisvaikutus RR-menetelmällä\n")
oprint(head(kokonaisvaikutusrr@output))
}
geombar <- function(object, title = "", y = "", x = "") {
out <- ggplot(object@output, aes(weight = Result, x = Tupakka, fill = Talotyyppi)) + geom_bar(position = "dodge") +
theme_grey(base_size = 24) +
labs( # label names
title = title,
y = y,
x = x
)
return(out)
}
ograph(RR, x = erottelu1, fill = erottelu2, title = 'Suhteellinen yksilöriski altistumattomaan verrattuna')
geombar(vaikutuspaf, title = "Vaikutus PAF-menetelmällä")
geombar(vaikutusrr, title = "Vaikutus RR-menetelmällä")
ggplot(väestö@output, aes(x = Talotyyppi, weight = Result)) + geom_bar() + labs(title = "Väestö")
ggplot(altistus@output[altistus@output$Altiste == "Radon" , ], aes(x = altistusResult, fill = Talotyyppi)) + geom_density(alpha = 0.2)
## POISTETTIIN KOKO KUNTAKARTTATOIMINNALLISUUS, KOSKA JOTKIN OSAT OVAT VANHENTUNEET. ON HELPOMPI RAKENTAA ALUSTA UUDESTAAN.
| |