PM25-HIA-methodology

example-calculations.R

library(units)
install_symbolic_unit("person")
install_conversion_constant("person", "death", const = -1)
options(digits = 8)

# For convenience, let's assume a population of 1 million people.
pop <- as_units(1e6, "person")

# Let's assume that the baseline annual all-cause mortality rate is 
# about 1%, i.e., about 10,000 per million (per year).
(y_all_0  <- as_units(10e3, "death/(Mperson*yr)")) 

# Let's assume that 1% of that is cause-specific (i.e., due to X).
(y_PM25_0  <- set_units(0.01 * y_all_0, "death/(Mperson*yr)"))

# Let's assume this is the cause-specific risk ratio, i.e., 
# the multiplicative change in risk per +10 ug/m3.
β <- 1.03

# Since our unit change in X is +10 ug/m3, this corresponds
# to an increase of +1 ug/m3, i.e., 0.1 * +10 ug/m3.
Δx <- 0.1

# As in Fang (2013) ... here is our first result.
show(y_all_0 * (1 - exp(-β * Δx)) * pop)

# As in Fann (2011) ... here is our second result. It should be the same as the first.
show(y_PM25_0 * (exp(β * Δx) - 1) * pop)