|
|
| (One intermediate revision by the same user not shown) |
| Line 3: |
Line 3: |
| | Tässä rakennetaan uudempaa versiota rakennuskantamallista. Lähtökohtana on kesäkuussa 2014 Rotterdamin Urgenche-kokoukseen kehitetty malli, johon tuli jokin bugi ja jolla ei lopulta saatu mitään tuloksia. Nyt (syyskuu 2014) bugit siivotaan ja mallin seuraava versio kehitetään. Lähtökohtana on [http://en.opasnet.org/en-opwiki/index.php?title=Climate_change_policies_and_health_in_Kuopio&oldid=33319#Results tämä malliversio] 25.8.2014. | | Tässä rakennetaan uudempaa versiota rakennuskantamallista. Lähtökohtana on kesäkuussa 2014 Rotterdamin Urgenche-kokoukseen kehitetty malli, johon tuli jokin bugi ja jolla ei lopulta saatu mitään tuloksia. Nyt (syyskuu 2014) bugit siivotaan ja mallin seuraava versio kehitetään. Lähtökohtana on [http://en.opasnet.org/en-opwiki/index.php?title=Climate_change_policies_and_health_in_Kuopio&oldid=33319#Results tämä malliversio] 25.8.2014. |
| | | | |
| − | <rcode graphics=1>
| + | Malli saatu taas pyörähtämään, ja se on palautettu alkuperäiselle paikalleen [[Climate change policies and health in Kuopio#Results]] 22.9.2014. Tällä paikalla olleen tilapäismallin uusin versio on [http://en.opasnet.org/en-opwiki/index.php?title=Talk:Climate_change_policies_and_health_in_Kuopio&oldid=33849#Ty.C3.B6versio_rakennusmallista]. |
| − | ############## TILAPÄISESTI HIA-FUNKTIOITA TÄHÄN
| |
| − | library(OpasnetUtils)
| |
| − | | |
| − | dose <- Ovariable("dose", # This calculates the body-weight-scaled exposure or "dose" to be used with ERFs.
| |
| − | dependencies = data.frame(Name = c(
| |
| − | "ERF", # Exposure-response function of the pollutants or agents (RR for unit exposure)
| |
| − | "exposure", # Exposure to the pollutants
| |
| − | "BW" # body weight
| |
| − | )),
| |
| − | formula = function(...) {
| |
| − | | |
| − | #################################################################33
| |
| − | ######### Body weight scaling: In some cases, exposure is given as per body weight and in some cases as absolute amounts.
| |
| − | # Here we add one index to define for this difference.
| |
| − |
| |
| − | scaling <- unique(ERF@output[c("ERF_parameter", "Exposure_agent")])
| |
| − |
| |
| − | scaling <- Ovariable("scaling", data = data.frame(
| |
| − | scaling,
| |
| − | Result = 1 * grepl("bw", scaling$ERF_parameter) # if ERF parameter contains bw, scaling is TRUE, i.e. 1.
| |
| − | ))
| |
| − | out <- exposure / (((BW - 1) * scaling) + 1) # If scaling is 0, BW cancels out.
| |
| − | out <- out * Ovariable( # Create alternative scenario with zero exposure.
| |
| − | output = data.frame(Exposcen = c("BAU", "No exposure"), Result = c(1, 0)),
| |
| − | marginal = c(TRUE, FALSE)
| |
| − | )
| |
| − |
| |
| − | return(out)
| |
| − | }
| |
| − | )
| |
| − | | |
| − | RR <- Ovariable("RR", # This calculates the total number of cases in each population subgroup.
| |
| − | # The cases are calculated for specific (combinations of) causes. However, these causes are NOT visible in the result.
| |
| − | dependencies = data.frame(Name = c(
| |
| − | "ERF", # Exposure-response function of the pollutants or agents (RR for unit exposure)
| |
| − | "dose", # Exposure to the pollutants
| |
| − | "frexposed", # fraction of population that is exposed
| |
| − | "bgexposure", # Background exposure (a level below which you cannot get in practice)
| |
| − | "threshold" # exposure level below which the agent has no impact.
| |
| − | )),
| |
| − | formula = function(...) {
| |
| − | | |
| − | ####################################################################
| |
| − | ####### This part is about risks relative to background.
| |
| − | # Calcualte the risk ratio to each subgroup based on the exposure in that subgroup.
| |
| − | # Combine pollutant-specific RRs by multiplying. For description, see [[Exposure-response function]].
| |
| − | test <- list()
| |
| − | marginals <- character()
| |
| − | | |
| − | #First take the relative risk estimates
| |
| − | ERFrr <- ERF
| |
| − | ERFrr@output <- ERF@output[ERF@output$ERF_parameter %in% c("RR", "RR bw") , ]
| |
| − | if(nrow((ERFrr * dose)@output) > 0) {
| |
| − | bigger <- (threshold > bgexposure) * threshold + (1 - (threshold > bgexposure)) * bgexposure # Choose bigger
| |
| − | RRrr <- exp(log(ERFrr) * (dose - bigger)) # Actual function
| |
| − | RRrr <- (RRrr - 1) * (dose > threshold) + 1 # RR is 1 below threshold
| |
| − | test <- c(test, RRrr)
| |
| − | marginals <- c(marginals, colnames(RRrr@output)[RRrr@marginal])
| |
| − | }
| |
| − |
| |
| − | # Then take the relative Hill estimates
| |
| − | ED50 <- threshold
| |
| − | ED50@output <- ED50@output[ED50@output$ERF_parameter %in% c("Relative Hill", "Relative Hill bw") , ]
| |
| − | if(nrow((ED50 * dose)@output) > 0 & nrow(ED50@output) > 0) {
| |
| − | RRrhill <- 1 + (dose * ERF) / (dose + ED50) # ERF has parameter value for Imax. If Imax < 0, risk reduces.
| |
| − | test <- c(test, RRrhill)
| |
| − | marginals <- c(marginals, colnames(RRrhill@output)[RRrhill@marginal])
| |
| − | }
| |
| − | | |
| − | # We need OR but not yet crucial, so let's postpone this. See [[:op_en:Converting between exposure-response parameters]]
| |
| − | # #Then take the odds ratio estimates
| |
| − | # OR <- ERF
| |
| − | # OR@output <- ERF@output[ERF@output$ERF_parameter %in% c("OR", "OR bw") , ]
| |
| − | # if(nrow((OR * dose)@output) > 0) {
| |
| − | # bigger <- (threshold > bgexposure) * threshold + (1 - (threshold > bgexposure)) * bgexposure # Choose bigger
| |
| − | # OR <- RR = OR/( 1-PX0+OR*PX0 ) # Actual function where PX0 is the background incidence. How to write a code?
| |
| − | # RRor <- (RRrr - 1) * (dose > threshold) + 1 # RR is 1 below threshold
| |
| − | # test <- c(test, RRor)
| |
| − | # marginals <- c(marginals, colnames(OR@output)[OR@marginal]
| |
| − | # }
| |
| − |
| |
| − | if(length(test) == 0) return(data.frame())
| |
| − | if(length(test) == 1) out <- test[[1]]@output
| |
| − | if(length(test) == 2) out <- orbind(test[[1]], test[[2]])
| |
| − | if(length(test) == 3) out <- orbind(orbind(test[[1]]), test[[2]], test[3])
| |
| − | | |
| − | out <- Ovariable(output = out, marginal = colnames(out) %in% marginals)
| |
| − | | |
| − | # Dilute the risk in the population if not all are exposed i.e. frexposed < 1.
| |
| − | out <- frexposed * (out - 1) + 1
| |
| − | out <- unkeep(out, prevresults = TRUE, sources = TRUE)
| |
| − | out <- oapply(out, cols = "Exposure_agent", FUN = prod)
| |
| − |
| |
| − | return(out)
| |
| − | }
| |
| − | )
| |
| − | | |
| − | totcases <- Ovariable("totcases", # This calculates the total number of cases in each population subgroup.
| |
| − | # The cases are calculated for specific (combinations of) causes. However, these causes are NOT visible in the result.
| |
| − | dependencies = data.frame(Name = c(
| |
| − | "population", # Population divided into subgroups as necessary
| |
| − | "disincidence", # Incidence of the disease of interest
| |
| − | "RR", # Relative risks for the given exposure
| |
| − | "ERF", # Other ERFs than those that are relative to background.
| |
| − | "dose", # Exposure to the pollutants
| |
| − | "frexposed", # fraction of population that is exposed
| |
| − | "threshold" # exposure level below which the agent has no impact.
| |
| − | )),
| |
| − | formula = function(...) {
| |
| − | | |
| − | test <- list()
| |
| − | marginals <- character()
| |
| − |
| |
| − | ############### First look at the relative risks based on RR
| |
| − | if(nrow((RR * dose)@output) > 0) {
| |
| − |
| |
| − | # if(class(population) != "ovariable") population <- EvalOutput(Ovariable("population", data = data.frame(Result = population)))
| |
| − | | |
| − | # # Remove redundant columns and locations.
| |
| − | # population@output <- dropall(population@output)
| |
| − | # population <- unkeep(population, sources = TRUE, prevresults = TRUE)
| |
| − | # disincidence <- unkeep(disincidence, sources = TRUE, prevresults = TRUE)
| |
| − |
| |
| − | # takeout is a vector of column names of indices that ARE in population but NOT in the disease incidence.
| |
| − | # However, populationSource is kept because oapply does not run if there are no indices.
| |
| − | | |
| − | if(class(population) == "ovariable") {
| |
| − | takeout <- setdiff(colnames(population@output), c(colnames(disincidence@output), "populationSource"))
| |
| − | pop <- oapply(population, cols = takeout, FUN = sum) # Aggregate to larger subgroups.
| |
| − | } else {
| |
| − | takeout <- character()
| |
| − | pop <- population
| |
| − | }
| |
| − | | |
| − | # pci is the proportion of cases across different population subroups based on differential risks and
| |
| − | # population sizes. pci sums up to 1 for each larger subgroup found in disincidence.
| |
| − | # See [[Population attributable fraction]].
| |
| − | pci <- population * RR
| |
| − | # Divide pci by the values of the actually exposed group (discard nonexposed)
| |
| − | temp <- pci
| |
| − | temp@output <- temp@output[temp@output$Exposcen == "BAU" , ]
| |
| − | temp <- unkeep(temp, cols = "Exposcen", prevresults = TRUE, sources = TRUE)
| |
| − | pci <- pci / oapply(temp, cols = takeout, FUN = sum)
| |
| − | | |
| − | # The cases are divided into smaller subgroups based on weights in pci.
| |
| − | # This is why the larger groups of population are used (pop instead of population).
| |
| − | out1 <- disincidence * pop * unkeep(pci, prevresults = TRUE, sources = TRUE)
| |
| − | out1 <- unkeep(out1, cols = "populationResult") # populationResult comes from pop and not from pci that actually contains
| |
| − | # the population weighting for takeout indices. Therefore it would be confusing to leave it there.
| |
| − | test <- c(test, out1)
| |
| − | marginals <- c(marginals, colnames(out1@output)[out1@marginal])
| |
| − | }
| |
| − | | |
| − | ##########################################################################
| |
| − | ############# This part is about absolute risks (i.e., risk is not affected by background rates).
| |
| − |
| |
| − | # Unit risk (UR), cancer slope factor (CSF), and Exposure-response slope (ERS) estimates.
| |
| − |
| |
| − | UR <- ERF
| |
| − | UR@output <- UR@output[UR@output$ERF_parameter %in% c("UR", "CSF", "ERS", "UR bw", "CSF bw", "ERS bw") , ]
| |
| − | if(nrow((UR * dose)@output) > 0) {
| |
| − | UR <- threshold + UR * dose * frexposed # Actual equation
| |
| − | # threshold is here interpreted as the baseline response (intercept of the line). It should be 0 for
| |
| − | # UR and CSF but it may have meaningful values with ERS
| |
| − | UR <- unkeep(UR, prevresults = TRUE, sources = TRUE)
| |
| − | UR <- oapply(UR, cols = "Exposure_agent", FUN = sum)
| |
| − | UR@output <- UR@output[!is.na(result(UR)) , ] # Remove empty rows
| |
| − | | |
| − | UR <- population * UR
| |
| − | test <- c(test, UR)
| |
| − | marginals <- c(marginals, colnames(UR)[UR@marginal])
| |
| − | }
| |
| − | | |
| − | # Step estimates: value is 1 below threshold and above ERF, and 0 in between.
| |
| − | # frexposed cannot be used with Step because this may be used at individual and maybe at population level.
| |
| − | Step <- ERF
| |
| − | Step@output <- Step@output[Step@output$ERF_parameter %in% c(
| |
| − | "Step", "Step bw", "ADI", "ADI bw", "TDI", "TDI bw", "RDI", "RDI bw", "NOAEL", "NOAEL bw") , ]
| |
| − | if(nrow((Step * dose)@output) > 0) {
| |
| − | Step <- 1 - (dose >= threshold) * (dose <= Step) # Actual equation
| |
| − | Step <- unkeep(Step, prevresults = TRUE, sources = TRUE)
| |
| − | | |
| − | Step <- oapply(Step, cols = "Exposure_agent", FUN = sum)
| |
| − | # out3 <- (population * 0 + 1) * Step # This is to maintain the ovariable structure # Does not work because population may be numeric.
| |
| − | test <- c(test, Step)
| |
| − | marginals <- c(marginals, colnames(Step@output)[Step@marginal])
| |
| − | }
| |
| − |
| |
| − | #####################################################################
| |
| − | # Combining effects
| |
| − | if(length(test) == 0) return(data.frame())
| |
| − | if(length(test) == 1) out <- test[[1]]
| |
| − | if(length(test) == 2) out <- orbind(test[[1]], test[[2]])
| |
| − | if(length(test) == 3) out <- orbind(orbind(test[[1]], test[[2]]), test[[3]])
| |
| − | | |
| − | out <- Ovariable(output = out, marginal = colnames(out) %in% marginals)
| |
| − |
| |
| − | return(out)
| |
| − | }
| |
| − | )
| |
| − | | |
| − | ################################## Varsinainen malli alkaa tästä
| |
| − | | |
| − | library(OpasnetUtils)
| |
| − | library(ggplot2)
| |
| − | library(rgdal)
| |
| − | library(maptools)
| |
| − | library(RColorBrewer)
| |
| − | library(classInt)
| |
| − | library(OpasnetUtilsExt)
| |
| − | library(RgoogleMaps)
| |
| − | | |
| − | language <- "EN"
| |
| − | openv.setN(1)
| |
| − | | |
| − | objects.latest("Op_en6007", code_name = "answer") # [[OpasnetUtils/Drafts]] findrest
| |
| − | objects.latest("Op_en6289", code_name = "initiate") # [[Building model]] buildings, heatingEnergy. # Generic building model.
| |
| − | objects.latest("Op_en5417", code_name = "initiate") # [[Population of Kuopio]] population
| |
| − | objects.latest("Op_en5932", code_name = "initiate") # [[Building stock in Kuopio]] Building ovariables: buildingStock, renovation,
| |
| − | #renovationShares, construction, constructionAreas, buildingTypes, heatingShares, heatingSharesNew, eventyear, efficiencies
| |
| − | objects.latest("Op_en3435", code_name = "disperse") # [[Exposure to PM2.5 in Finland]] iF, emissionLocations
| |
| − | #Summarised Piltti matrix, another copy of the code on a more reasonable page
| |
| − | # Default run: en.opasnet.org/en-opwiki/index.php?title=Special:RTools&id=aXDIVDboftr1bTEd
| |
| − | objects.latest("Op_en2791", code_name = "initiate") # [[Emission factors for burning processes]] emissionFactors, fuelTypes
| |
| − | objects.latest("Op_en5488", code_name = "initiate") # [[Energy use of buildings]] energyUse, efficienciesNew, savingPotential,
| |
| − | #objects.latest('Op_en2261', code_name = 'initiate') # [[Health impact assessment]] ovariables totcases, AF, totcases2.
| |
| − | objects.latest('Op_en5917', code_name = 'initiate') # [[Disease risk]] disincidence
| |
| − | objects.latest('Op_en5827', code_name = 'initiate') # [[ERFs of environmental pollutants]] ERF, threshold
| |
| − | #objects.latest('Op_en5453', code_name = 'initiate') # [[Burden of disease in Finland]] BoD
| |
| − | | |
| − | frexposed <- 1 # fraction of population that is exposed
| |
| − | bgexposure <- 0 # Background exposure to an agent (a level below which you cannot get in practice)
| |
| − | BW <- 70 # Body weight (is needed for RR calculations although it is irrelevant for PM2.5)
| |
| − | | |
| − | directs <- tidy(opbase.data("Op_en5461", subset = "Direct inputs"), direction = "wide")
| |
| − | colnames(directs) <- gsub(" ", "_", colnames(directs))
| |
| − | | |
| − | colnames(iF@output)[colnames(iF@output) == "City.area"] <- "Emission_site"
| |
| − | colnames(iF@output)[colnames(iF@output) == "Emissionheight"] <- "Emission_height"
| |
| − | colnames(emissionLocations@data) <- gsub(" ", "_", colnames(emissionLocations@data))
| |
| − | | |
| − | # Calculate intake fractions for Kuopio districts
| |
| − | districts <- tidy(opbase.data("Op_en3435.kuopio_city_districts"), widecol = "Location") # [[Exposure to PM2.5 in Finland]]
| |
| − | colnames(districts) <- gsub("\\.", "_", colnames(districts))
| |
| − | districts <- Ovariable("districts", data = data.frame(districts, Result = 1))
| |
| − | | |
| − | if(language == "EN") deci <- "Decisions" else deci <- "Päätökset"
| |
| − | | |
| − | decisions <- opbase.data('Op_en5461', subset = deci) # [[Climate change policies and health in Kuopio]]
| |
| − | | |
| − | DecisionTableParser(decisions)
| |
| − | | |
| − | # Remove previous decisions, if any.
| |
| − | rm(
| |
| − | "buildings",
| |
| − | "buildingStock",
| |
| − | "buildingTypes",
| |
| − | "construction",
| |
| − | "efficiencies",
| |
| − | "efficienciesNew",
| |
| − | "energyUse",
| |
| − | "heatingShares",
| |
| − | "heatingSharesNew",
| |
| − | "renovation",
| |
| − | "renovationShares",
| |
| − | "fuelTypes",
| |
| − | "year",
| |
| − | envir = openv
| |
| − | )
| |
| − | | |
| − | fuelTypes <- CheckDecisions(EvalOutput(fuelTypes, verbose = TRUE))
| |
| − | | |
| − | #################### Manage the data before calculating
| |
| − | | |
| − | # The building stock is measured as m^2 floor area.
| |
| − | | |
| − | if(!exists("language")) language <- "FI"
| |
| − | | |
| − | obsyear <- (192:205) * 10
| |
| − | | |
| − | ####### Remove columns and rows that are not needed from data. These must be done before EvalOutput.
| |
| − | | |
| − | # This is the population of Kuopio (i.e. population living in the building stock)
| |
| − | | |
| − | population <- EvalOutput(population, verbose = TRUE)
| |
| − | pop <- oapply(population, cols = "City_area", FUN = "sum") # Sum across city areas.
| |
| − | pop <- unkeep(pop * 1, prevresults = TRUE) # Remove populationResult column because it would cause trouble in merge.
| |
| − | pop <- population / pop
| |
| − | | |
| − | buildingStock@data <- buildingStock@data[
| |
| − | buildingStock@data$Observation == "AreaHR" ,
| |
| − | colnames(buildingStock@data) != "Observation"
| |
| − | ]
| |
| − | buildingStock <- buildingStock * pop
| |
| − | | |
| − | construction@data <- construction@data[
| |
| − | construction@data$Observation == "Area" ,
| |
| − | colnames(construction@data) != "Observation"
| |
| − | ]
| |
| − | construction <- construction * constructionAreas / 100 / 3 # Statistics are for three years (2010-2012)
| |
| − | | |
| − | energyUse@data <- energyUse@data[
| |
| − | energyUse@data[["Energy use"]] == "Heat" ,
| |
| − | colnames(energyUse@data) != "Energy use"
| |
| − | ]
| |
| − | | |
| − | savingPotential@data <- savingPotential@data[
| |
| − | savingPotential@data$Observation == "Relative" ,
| |
| − | colnames(savingPotential@data) != "Observation"
| |
| − | ]
| |
| − | savingPotential <- 1 - savingPotential / 100
| |
| − | # Note: If you add decisions to savingPotential, this formula must be updated.
| |
| − | | |
| − | # The data are not yet specific to construction year, so remove index:
| |
| − | heatingSharesNew@data <- heatingSharesNew@data[colnames(heatingSharesNew@data) != "Constructed"]
| |
| − | | |
| − | # Fill in Heating types and convert from % to fraction.
| |
| − | heatingSharesNew <- findrest(heatingSharesNew, cols = "Heating", total = 100)
| |
| − | result(heatingSharesNew) <- result(heatingSharesNew) / 100
| |
| − | | |
| − | heatingShares <- EvalOutput(heatingShares, verbose = TRUE)
| |
| − | result(heatingShares) <- result(heatingShares) / 100 # From % to fraction
| |
| − | | |
| − | # Fill in the rest of the data (the last emission class was omitted because it is determined by the total).
| |
| − | efficienciesNew <- findrest(efficienciesNew, cols = "Efficiency", total = 100)
| |
| − | result(efficienciesNew) <- result(efficienciesNew) / 100
| |
| − | | |
| − | # When renovations are done, which type are they?
| |
| − | renovationShares <- findrest(renovationShares, cols = "Renovation", total = 100)
| |
| − | result(renovationShares) <- result(renovationShares) / 100
| |
| − | | |
| − | # What fraction of buildings is renovated per year?
| |
| − | renovation <- EvalOutput(renovation, verbose = TRUE)
| |
| − | result(renovation) <- result(renovation) / 100 # From % to fraction
| |
| − | | |
| − | disincidence <- disincidence / 100000 # Change units from 1/100000py to 1/py.
| |
| − | | |
| − | ###################### Actual model
| |
| − | | |
| − | buildings <- EvalOutput(buildings, verbose = TRUE)
| |
| − | | |
| − | ####### Sort factors
| |
| − | | |
| − | buildings@output <- buildings@output[
| |
| − | buildings@output$EfficiencyPolicy == "BAU" ,
| |
| − | colnames(buildings@output) != "EfficiencyPolicy"
| |
| − | ]
| |
| − | | |
| − | buildings@output$Renovation <- factor(
| |
| − | buildings@output$Renovation,
| |
| − | levels = c("None", "General", "Windows", "Technical systems", "Sheath reform"),
| |
| − | ordered = TRUE
| |
| − | )
| |
| − | | |
| − | buildings@output$Efficiency <- factor(
| |
| − | buildings@output$Efficiency,
| |
| − | levels = c("Old", "New", "Low-energy", "Passive"),
| |
| − | ordered = TRUE
| |
| − | )
| |
| − | | |
| − | buildings@output$RenovationPolicy <- factor(
| |
| − | buildings@output$RenovationPolicy,
| |
| − | levels = c("BAU", "Active renovation", "Effective renovation"),
| |
| − | ordered = TRUE
| |
| − | )
| |
| − | | |
| − | heatingEnergy <- EvalOutput(heatingEnergy, verbose = TRUE)
| |
| − | heatingEnergy <- unkeep(heatingEnergy, cols = "Building2", sources = TRUE, prevresults = TRUE)
| |
| − | | |
| − | if(exists("server")) # The following code only works from Opasnet server.
| |
| − | {
| |
| − | temp <- buildings * districts
| |
| − | temp@output <- temp@output[temp@output$Year == 2040 , ]
| |
| − | temp <- unkeep(temp, sources = TRUE, prevresults = TRUE)
| |
| − | temp@output <- dropall(temp@output)
| |
| − | temp <- oapply(temp, cols = c("Building", "Heating", "Efficiency", "Renovation"), FUN = "sum", na.rm = TRUE)
| |
| − | | |
| − | if(language == "EN") lege <- "Floor area" else lege <- "Rakennuspinta-ala"
| |
| − | | |
| − | MyRmap(
| |
| − | ova2spat(temp, coord = c("E", "N"), proj4string = "+init=epsg:3067"), # National Land Survey uses EPSG:3067 (ETRS-TM35FIN)
| |
| − | plotvar = "Result", legend_title = lege, numbins = 8, pch = 19, cex = 2
| |
| − | )
| |
| − | }
| |
| − | | |
| − | ####### Calculate emissions
| |
| − | | |
| − | emis <- heatingEnergy
| |
| − | emis@output <- emis@output[emis@output$Year >= 1980 , ]
| |
| − | emis <- oapply(emis, cols = c("Efficiency", "Renovation", "Building"), FUN = "sum")
| |
| − | emis@output$Year <- as.numeric(as.character(emis@output$Year))
| |
| − | | |
| − | emis <- emis * fuelTypes * emissionFactors * 3.6 * 1E-9 # convert from kWh /a to MJ /a and mg to ton
| |
| − | | |
| − | emis <- unkeep(emis, sources = TRUE, prevresults = TRUE)
| |
| − | | |
| − | emis <- emis * emissionLocations
| |
| − | emis@output$Emission_site <- ifelse(
| |
| − | emis@output$Emission_site == "At site of consumption",
| |
| − | as.character(emis@output$City_area),
| |
| − | as.character(emis@output$Emission_site)
| |
| − | )
| |
| − | | |
| − | emis2 <- oapply(
| |
| − | emis,
| |
| − | cols = c("Heating", "Burner", "Fuel", "City_area", "Emission_site", "Emission_height", "emissionLocationsSource"),
| |
| − | FUN = sum
| |
| − | )
| |
| − | | |
| − | ggplot(emis@output, aes(x = Year, weight = Result, fill = Emission_site)) + geom_bar() + facet_grid( Pollutant ~ ., scales = "free_y")
| |
| − | | |
| − | ### Use these population and iF values in health impact assessment. Why?
| |
| − | | |
| − | pop <- 5E+5
| |
| − | | |
| − | iF <- 1E-7
| |
| − | | |
| − | exposure <- Ovariable("exposure",
| |
| − | # emis is in ton /a
| |
| − | # iF = conc (g /m3) * pop (#) * BR (m3 /s) / emis (g /s) <=> conc = emis * iF / BR / pop # conc is the exposure concentration
| |
| − | dependencies = data.frame(Name = c("emis2", "iF", "pop")),
| |
| − | formula = function(...) {
| |
| − | BR <- 20 # Nominal breathing rate (m^3 /d)
| |
| − | BR <- BR / 24 / 3600 # m^3 /s
| |
| − | out <- 1E+12 / 365 / 24 / 3600 # Emission scaling from ton /a to ug /s.
| |
| − | out <- (emis2 * out) * iF / BR / pop # the actual equation
| |
| − | out <- unkeep(out, prevresults = TRUE, sources = TRUE)
| |
| − | return(out)
| |
| − | }
| |
| − | )
| |
| − | exposure <- EvalOutput(exposure, verbose = TRUE)
| |
| − | | |
| − | colnames(exposure@output)[colnames(exposure@output) == "Pollutant"] <- "Exposure_agent"
| |
| − | ggplot(exposure@output, aes(x = Year, weight = exposureResult)) + geom_bar() + facet_grid(Exposure_agent ~ FuelPolicy, scales = "free_y")
| |
| − | | |
| − | totcases <- EvalOutput(totcases)
| |
| − | #totcases@output <- totcases@output[!is.na(result(totcases)) , ] # Not needed with new oapply
| |
| − | | |
| − | cases <- totcases * Ovariable("incre", data = data.frame(Exposcen = c("BAU", "No exposure"), Result = c(1, -1)))
| |
| − | cases <- oapply(cases, cols = c("Exposcen", "Age", "Sex", "City_area"), FUN = sum)
| |
| − | | |
| − | DW <- Ovariable("DW", data = data.frame(directs["Trait"], Result = directs$DW))
| |
| − | L <- Ovariable("L", data = data.frame(directs["Trait"], Result = directs$L))
| |
| − | | |
| − | #ggplot(cases@output[cases@output$Causes == "PM2.5+" , ], aes(x = Year, weight = Result/10, fill = Disease)) + geom_bar() + facet_grid(FuelPolicy ~ RenovationPolicy, scale = "free_y")
| |
| − | | |
| − | DALYs <- cases * DW * L
| |
| − | #DALYs <- oapply(DALYs, cols = c("Trait"), FUN = sum)
| |
| − | DALYs@output <- DALYs@output[DALYs@output$Trait != "Lung cancer" , ] # Has to be removed to avoid double counting.
| |
| − | | |
| − | #ggplot(DALYs@output[DALYs@output$Causes == "PM2.5+" , ], aes(x = Year, weight = Result/10, fill = Response.metric)) + geom_bar() + facet_grid(FuelPolicy ~ RenovationPolicy, scale = "free_y")
| |
| − | ggplot(DALYs@output, aes(x = Year, weight = Result, fill = Trait)) + geom_bar() + facet_grid(RenovationPolicy ~ FuelPolicy)
| |
| − | | |
| − | | |
| − | ############################# Plot graphs
| |
| − | | |
| − | ######## Buildings in Kuopio on map
| |
| − | | |
| − | if(exists("server")) # The following code only works from Opasnet server.
| |
| − | {
| |
| − | if(language == "EN") lege2 <- "Age of building" else lege2 <- "Rakennuksen ikä"
| |
| − | | |
| − | shp <- readOGR('PG:host=localhost user=postgres dbname=spatial_db','kuopio_house')
| |
| − | proj4string(shp) <- ("+init=epsg:3067") # The map projection of NLS Finland.
| |
| − | | |
| − | epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
| |
| − | shp2 <- spTransform(shp,epsg4326String) # Convert to longitude-latitude projection.
| |
| − | | |
| − | MyRmap(shp2, plotvar = "ika", legend_title = lege2, numbins = 8, pch = 19, cex = 0.3) # Draw the map.
| |
| − | }
| |
| − | | |
| − | if(exists("server")) BS <- 24 else BS <- 12
| |
| − | | |
| − | ### Impact graphs in ENGLISH
| |
| − | | |
| − | if(language == "EN") {
| |
| − | | |
| − | plo <- ggplot(buildings@output) + geom_bar() + facet_grid(. ~ RenovationPolicy) + theme_gray(base_size = BS) +
| |
| − | aes(x = Year, weight = buildingsResult/1000000, fill = Renovation) + labs(y = "Floor area (M m2)",
| |
| − | title = "Building impacts of renovation policy")
| |
| − | | |
| − | print(plo)
| |
| − | | |
| − | plo <- ggplot(heatingEnergy@output) + geom_bar() + facet_grid(. ~ RenovationPolicy) + theme_gray(base_size = BS) +
| |
| − | labs(title = "Energy impacts of renovation policy", y = "Heating energy need (GWh /a)") +
| |
| − | aes(x = Year, weight = heatingEnergyResult/1E+6, fill = Heating)
| |
| − | | |
| − | print(plo)
| |
| − | | |
| − | emis@output <- emis@output[emis@output$Renovation == "BAU" , ]
| |
| − | | |
| − | plo <- ggplot(emis@output) + geom_bar() + facet_grid(Pollutant ~ FuelPolicy, scales = "free_y") + theme_gray(base_size = BS) +
| |
| − | labs(title = "Emission impacts of biofuel policy", y = "Emissions to air (ton /a)") +
| |
| − | aes(x = Year, weight = Result, fill = Fuel)
| |
| − | | |
| − | print(plo)
| |
| − | | |
| − | plo <- ggplot(health@output) + geom_bar() + facet_grid(FuelPolicy ~ RenovationPolicy) + theme_gray(base_size = BS) +
| |
| − | labs(title = "Health impacts of fuel and renovation policy", y = "Premature deaths (# /a)") +
| |
| − | aes(x = Year, weight = Result, fill = Heating)
| |
| − | | |
| − | print(plo)
| |
| − | }
| |
| − | | |
| − | ### Impact graphs in FINNISH
| |
| − | | |
| − | if(language == "FI") {
| |
| − | | |
| − | plo <- ggplot(buildings@output) + geom_bar() + facet_grid(. ~ Remonttipolitiikka) + theme_gray(base_size = BS) +
| |
| − | aes(x = Year, weight = buildingsResult/1000000, fill = Renovation) + labs(y = "Rakennuspinta-ala (M m2)",
| |
| − | title = "Remonttipolitiikan vaikutus rakennuskantaan")
| |
| − | | |
| − | print(plo)
| |
| − | | |
| − | plo <- ggplot(heatingEnergy@output) + geom_bar() + facet_grid(. ~ Remonttipolitiikka) + theme_gray(base_size = BS) +
| |
| − | labs(title = "Remonttipolitiikan vaikutus energiankulutukseen", y = "Lämmitysenergian tarve (GWh /a)") +
| |
| − | aes(x = Year, weight = heatingEnergyResult/1E+6, fill = Heating)
| |
| − | | |
| − | print(plo)
| |
| − | | |
| − | emis@output <- emis@output[emis@output$Remonttipolitiikka == "BAU" , ]
| |
| − | | |
| − | plo <- ggplot(emis@output) + geom_bar() + facet_grid(Pollutant ~ Bioenergiapolitiikka, scales = "free_y") +
| |
| − | theme_gray(base_size = BS) + labs(title = "Bioenergiapolitiikan päästövaikutukset", y = "Päästöt ilmaan (ton /a)") +
| |
| − | aes(x = Year, weight = Result, fill = Fuel)
| |
| − | | |
| − | print(plo)
| |
| − | | |
| − | plo <- ggplot(health@output) + geom_bar() + facet_grid(Bioenergiapolitiikka ~ Remonttipolitiikka) + theme_gray(base_size = BS) +
| |
| − | labs(title = "Remontti- ja bioenergiapolitiikan terveysvaikutukset", y = "Ennenaikaisia kuolemia (# /a)") +
| |
| − | aes(x = Year, weight = Result, fill = Heating)
| |
| − | | |
| − | print(plo)
| |
| − | }
| |
| − | | |
| − | ggplot(dose@output, aes(x = Year, y = doseResult, colour = RenovationPolicy)) + geom_point() + facet_grid(Exposcen ~ FuelPolicy)
| |
| − | | |
| − | ggplot(totcases@output, aes(x = Year, y = totcasesResult, colour = Exposcen)) + geom_point() + facet_grid(Trait ~ FuelPolicy)
| |
| − | | |
| − | ggplot(RR@output, aes(x = RRResult, colour = Trait)) + geom_density()
| |
| − | | |
| − | ggplot(cases@output, aes(x = Year, weight = Result, fill = Trait)) + geom_bar() + facet_grid(RenovationPolicy ~ FuelPolicy)
| |
| − | | |
| − | </rcode>
| |
| | | | |
| | ==Työtehtäviä Urgenche-kuukausikokouksesta== | | ==Työtehtäviä Urgenche-kuukausikokouksesta== |
Tässä rakennetaan uudempaa versiota rakennuskantamallista. Lähtökohtana on kesäkuussa 2014 Rotterdamin Urgenche-kokoukseen kehitetty malli, johon tuli jokin bugi ja jolla ei lopulta saatu mitään tuloksia. Nyt (syyskuu 2014) bugit siivotaan ja mallin seuraava versio kehitetään. Lähtökohtana on tämä malliversio 25.8.2014.
Kokouksessa esiteltiin URGENCHE- hanke ja sen taustat ja tavoitteet. Tämän jälkeen keskusteltiin energiantuotannon päästövähennysmahdollisuuksista Kuopiossa sekä päätettiin mitä asioita energiantuotantoa koskevaan päästövähennysskenaarioon otetaan mukaan. Päädyttiin selvittämään ja keräämään tietoa kolmesta asiakokonaisuudesta, jotka ovat:
Peter Seppälä (Kuopion Energia) arvioi, että Haapaniemi 3- ja Haapaniemi 2- voimalaitosten polttoaineina käytetään vuodesta 2014 eteenpäin 50 % turvetta ja 50 % biopolttoaineita. Tämä on polttoaineiden saatavuuden ja huoltovarmuuden näkökulmasta realistinen arvio, kun Haapaniemi 2- voimalaitos on muutettu pölypoltosta leijupoltolle. Tällöin Haapaniemen voimalaitosten CO2-päästöt olisivat noin 324 000 t vuodessa. Kyseisellä päästötasolla Kuopion ilmastopoliittisen ohjelman – 40 % päästötavoite vuoteen 1990 verrattuna on mahdollista saavuttaa jo vuonna 2014.
Sovittiin, että Marjo Niittynen laatii taulukon energiantuotantoskenaarioon tarvittavista tiedoista (vuodet, energia, polttoaineet, päästöt). Kuopion Energia täyttää taulukossa kysytyt tiedot tämänhetkisten arvioiden mukaan. Taulukossa on arviot mm. keskimääräisistä odotettavista polttoaineosuuksista ja -määristä, teoreettinen maksimi biopolttoaineille – vaihtoehto sekä näitä koskevat päästötiedot. Kuopion Energia toimittaa URGENCHE-hankkeelle arvion/mallin siitä, miten muutokset energiankulutuksessa (Haapaniemellä tuotettu sähkö/lämpö) tulee huomioida koko kaupungin CO2-päästöjä laskettaessa.
THL:lla keskustellaan vielä laaditaanko tällainen skenaario vai ei. Ajatuksena olisi laskea, millainen muutos Haapaniemen energiantuotannossa (ja päästöissä) tapahtuisi, jos XX % Kuopion pien- ja kerrostalojen lämmityksestä siirtyisi kiinteistökohtaiseen lämmitykseen biopolttoaineilla. Lisäksi arvioidaan millainen vaikutus tällä olisi koko kaupungin päästöihin (CO2 ja hiukkaset). Jos päästötiedot aiotaan saada myös paikkatietona (leviämismalli), tästä pitää tehdä nopea päätös, että tulokset saadaan Thessalonikin kokoukseen huhtikuussa.
Hajautettua energiantuotantoa ei tutkita aluelämmön (esim. pienen alueen hakekattilat) käytön näkökulmasta sen epärealistisuuden vuoksi.
Selvitetään, miten laajamittainen (=jokin arvio) siirtyminen maalämmön, aurinkoenergian tai tuulivoiman käyttöön vaikuttaisi Haapaniemen lämpökuorman tarpeeseen ja päästöihin sekä koko kaupungin päästöihin. Asiaa selvitetään ja tehdään laskennallinen arvio.
Kokouksessa esiteltiin URGENCHE- hanke ja sen taustat ja tavoitteet. Tämän jälkeen keskusteltiin liikenteen ja maankäytön päästövähennysmahdollisuuksista Kuopiossa sekä päätettiin mitä asioita liikennettä ja maankäyttöä koskevaan päästövähennysskenaarioon otetaan mukaan.
Maankäyttöön ja liikenteeseen liittyen on eri tahoilla ajatuksia ja suunnitelmia, jotka eivät välttämättä ole lähellä toteutumista, mutta jotka voivat olla URGENCHEN kannalta mielenkiintoisia. Näiden osalta voidaan URGENCHE-hankkeessa tehdä arviointeja, jotka antavat tietoa toteuttamisen seurauksista tai toteutusvaihtoehtojen eroista.
Joukkoliikenne kilpailutetaan huomattavasti nykyistä laajemmin ensi vuoden aikana. Kuopiossa joukkoliikenteen käyttäjämäärissä on ollut 4-5 % kasvua tänä vuonna.
Kokouksessa esiteltiin URGENCHE- hanke ja sen taustat ja tavoitteet. Tämän jälkeen keskusteltiin rakennusten energiatehokkuuteen liittyvistä päästövähennysmahdollisuuksista Kuopiossa sekä päätettiin miten edetään tätä asiakokonaisuutta koskevan päästövähennysskenaarion kanssa.
