Assessment of the health impacts of H1N1 vaccination: Difference between revisions

From Opasnet
Jump to navigation Jump to search
Line 117: Line 117:
# however if the probability has an uncertainty, only 1 sample from each distribution is sampled where
# however if the probability has an uncertainty, only 1 sample from each distribution is sampled where
# samples of the probability will be used as inputs; disabled sampling here since it would make different scenarios incomparable
# samples of the probability will be used as inputs; disabled sampling here since it would make different scenarios incomparable
tlcf <- data.frame(Result = 0.1)#obs = 1:n, Result=rbeta(n, paralpha(0.1, 0.05^2), parbeta(0.1, 0.05^2)))


sf <- function(isfp = sfp(), ipop = op_baseGetData("opasnet_base", "Op_en2949")[,-c(1,2)], iimm = imm(), n = 1000,  
sf <- function(isfp = sfp(), ipop = op_baseGetData("opasnet_base", "Op_en2949")[,-c(1,2)], iimm = imm(), n = 1000,  
Line 149: Line 147:
# distribution. But we're interested in the uncertainty of the estimate... When n is suitably large a poisson distribution becomes  
# distribution. But we're interested in the uncertainty of the estimate... When n is suitably large a poisson distribution becomes  
# an excellent approximation of the binomial distribution. So we can fit the poisson distribution using an R function (I didn't find
# an excellent approximation of the binomial distribution. So we can fit the poisson distribution using an R function (I didn't find
# one for binomial) from just one observation -> acquire value of x with standard error -> use as uncertainty of estimate.
# one for binomial) from just one observation -> acquire value of x with standard error -> use as uncertainty of estimate; which produces
# rubbish because we only have a single observation, sd of the mean estimate of the poisson distribution appears to be (n*p)^0.5, so the
# variance of the parameter is n * p, which divided by n gives p ... so the expected value is p and the variance is p, which causes problems
# with the beta distribution (shape parameters become negative), hence variance is assumed to be variance^2 arbitrarily


sfp <- function(dsf = op_baseGetData("opasnet_base", "Op_en4933")[,-c(1,2)], dpop = op_baseGetData("opasnet_base", "Op_en2949",  
sfp <- function(dsf = op_baseGetData("opasnet_base", "Op_en4933")[,-c(1,2)], dpop = op_baseGetData("opasnet_base", "Op_en2949",  
Line 164: Line 165:
for (i in 1:nrow(out)) {
for (i in 1:nrow(out)) {
temp2 <- data.frame(obs=1:n, Age=out$Age[i], fitp=rbeta(n, paralpha(fitdistr(round(out$dsf[i]), "poisson")[[1]] /  
temp2 <- data.frame(obs=1:n, Age=out$Age[i], fitp=rbeta(n, paralpha(fitdistr(round(out$dsf[i]), "poisson")[[1]] /  
out$popnonimm[i], fitdistr(round(out$dsf[i]), "poisson")[[2]]^2 / out$popnonimm[i]^2),  
out$popnonimm[i], (fitdistr(round(out$dsf[i]), "poisson")[[2]]^2 / out$popnonimm[i])^2),  
parbeta(fitdistr(round(out$dsf[i]), "poisson")[[1]] / out$popnonimm[i], fitdistr(round(out$dsf[i]),  
parbeta(fitdistr(round(out$dsf[i]), "poisson")[[1]] / out$popnonimm[i], (fitdistr(round(out$dsf[i]),  
"poisson")[[2]]^2 / out$popnonimm[i]^2)))
"poisson")[[2]]^2 / out$popnonimm[i])^2)))
out2 <- rbind(out2, merge(out[i,], temp2, all.x = TRUE, all.y = FALSE))
out2 <- rbind(out2, merge(out[i,], temp2, all.x = TRUE, all.y = FALSE))
}
}
Line 191: Line 192:
out <- model.frame(I(dfrg * dsf / ilcf) ~., data = out)
out <- model.frame(I(dfrg * dsf / ilcf) ~., data = out)
colnames(out)[1] <- "Result"
colnames(out)[1] <- "Result"
out <- data.frame(Result=dmsf$dmsf/sum(out$Result))
if (sum(as.numeric(colnames(out)=="obs"))>0) {#out <- as.data.frame(as.table(tapply(out$)))
out <- data.frame(obs = levels(factor(out$obs)), Result = dmsf$dmsf/tapply(out$Result, out$obs, sum))
} else {out <- data.frame(Result=dmsf$dmsf/sum(out$Freq))}
return(out)
return(out)
}
}
Line 215: Line 218:
library(MASS)
library(MASS)
library(ggplot2)
library(ggplot2)
n <- 1000
tlcf <- data.frame(obs = 1:n, Result=rbeta(n, paralpha(0.2, 0.1^2), parbeta(0.2, 0.1^2)))
#ggplot(tlcf, aes(x=Result)) + geom_density()
#tlcf <- data.frame(Result = 0.1)
#test2 <- data.frame(obs = 1:n, Result=rbeta(n, paralpha(0.2, 0.2^2+0.01), parbeta(0.2, 0.2^2+0.01)))
#ggplot(test2, aes(x=Result)) + geom_density()


# Default, model based on current knowledge
# Default, model based on current knowledge
Line 237: Line 248:
temp <- as.data.frame(as.table(tapply(test$Result, test[,c("obs","Scenario")], sum)))
temp <- as.data.frame(as.table(tapply(test$Result, test[,c("obs","Scenario")], sum)))
ggplot(temp, aes(x = Freq, fill = Scenario)) + geom_density(alpha = 0.2)
ggplot(temp, aes(x = Freq, fill = Scenario)) + geom_density(alpha = 0.2)
tapply(temp$Freq, temp[,c("Scenario")], mean)


#vacscenario <- outcome()
#vacscenario <- outcome()

Revision as of 13:55, 5 April 2011



Scope

What was the overall health impact of the H1N1 vaccination, given that there is some evidence of it having caused narcolepsy.

Definition

Causal diagram.
Variables
Indicators

Formula

  • Basic model, works.
    • Uncertainties of ERF of vaccine on narcolepsy and probability of catching swine flu are implemented.

+ Show code

Result

  • From initial results it would appear like swine flu is more significant than narcolepsy in terms of DALYs.
    • Vaccinating as planned would result in approximately 2500 DALYs due to swine flu and narcolepsy combined.
    • Vaccinating no-one would result in approximately 4000 DALYs due to swine flu.

Results

Conclusions

See also

References

  1. ↑ WHO life table estimates for 2008, Finland
Retrieved from "https://en.opasnet.org/index.php?title=Assessment_of_the_health_impacts_of_H1N1_vaccination&oldid=19971"