diff --git a/R/plotValuationTables.R b/R/plotValuationTables.R index e056851d5dc0718de1f09751ed6917f15871e80a..1e05048be225f4c43de801e4d34331d9ce98c038 100644 --- a/R/plotValuationTables.R +++ b/R/plotValuationTables.R @@ -1,20 +1,67 @@ #' @export -makeQxDataFrame = function(..., YOB=1972, Period=NA) { - data=list(...); - names(data) = lapply(data, function(t) t@name); - if (missing(Period)) { - cat("Year of birth: ", YOB, "\n"); - data = lapply(data, function(t) cbind(x=ages(t), y=deathProbabilities(t, YOB=YOB))) - } else { - cat("Period: ", Period,"\n"); - data = lapply(data, function(t) cbind(x=ages(t), y=periodDeathProbabilities(t, Period=Period))) - } +#' If reference is given, normalize all probabilities by that table! +makeQxDataFrame = function(..., YOB=1972, Period=NA, reference=NULL) { + data=list(...); + names(data) = lapply(data, function(t) t@name); + reference_ages = NULL; + + if (missing(Period)) { + cat("Year of birth: ", YOB, "\n"); + if (!missing(reference)) { + reference_ages = ages(reference); + reference = deathProbabilities(reference, YOB=YOB); + } + data = lapply(data, function(t) { + str(deathProbabilities(t, YOB=YOB)); + normalize_deathProbabilities( + cbind(x=ages(t), y=deathProbabilities(t, YOB=YOB)), + reference = reference, + referenceAges = reference_ages) + }); + } else { + cat("Period: ", Period,"\n"); + if (!missing(reference)) { + reference_ages = ages(reference); + reference = periodDeathProbabilities(reference, Period=Period); + } + data = lapply(data, function(t) { + str(deathProbabilities(t, YOB=YOB)); + normalize_deathProbabilities( + cbind(x=ages(t), y=periodDeathProbabilities(t, Period=Period)), + reference = reference, + referenceAges = reference_ages) + }); + } + + list.names = names(data) + lns <- sapply(data, nrow) + data <- as.data.frame(do.call("rbind", data)) + data$group <- rep(list.names, lns) + data +} + +normalize_deathProbabilities = function (data, reference=NULL, referenceAges=NULL) { + if (missing(reference) || missing(referenceAges) || is.null(reference) || is.null(referenceAges)) { + return(data); + } +# browser(); + # Find which ages exist in both and obtain those indices from the data and the reference list: + useages = intersect(data[,"x"], referenceAges) + dataindices = match(useages, data[,"x"]) + refindices = match(useages, referenceAges) + + # Find which ages in data do NOT exist in the reference ages (and are thus NOT normalized at all) + # Print a warning! + missingrefs = setdiff(data[,"x"], referenceAges) + if (length(missingrefs)>0) { + warning("Reference mortality table does not contain ages ", + missingrefs, + " required for normalization. These ages will not be normalized!") + } - list.names = names(data) - lns <- sapply(data, nrow) - data <- as.data.frame(do.call("rbind", data)) - data$group <- rep(list.names, lns) - data + # Now divide the data by the corresponding entries from the reference list + data[dataindices, "y"] = data[dataindices, "y"] / reference[refindices] + data } #' Plot multiple valuation tables (life tables) in one plot @@ -28,7 +75,7 @@ makeQxDataFrame = function(..., YOB=1972, Period=NA) { #' @param legend.key.width The keywith of the lines in the legend (default is \code{unit(25,"mm")}) #' #' @export -plotValuationTables = function(data, ..., title = "", legend.position=c(0.9,0.1), legend.key.width = unit(25, "mm")) { +plotValuationTables = function(data, ..., xlim=NULL, ylim=NULL, title = "", legend.position=c(0.9,0.1), legend.key.width = unit(25, "mm")) { if (!is.data.frame(data)) { data = makeQxDataFrame(data, ...); } @@ -60,6 +107,7 @@ plotValuationTables = function(data, ..., title = "", legend.position=c(0.9,0.1) #labels = scales::trans_format('log10', scales::math_format(10^.x)) ) + + coord_cartesian(xlim=xlim, ylim=ylim) + annotation_logticks(sides="lr") + xlab("Alter") + labs(colour="Sterbetafel"); if (title != "") { @@ -67,6 +115,62 @@ plotValuationTables = function(data, ..., title = "", legend.position=c(0.9,0.1) } pl } + +plotValuationTableComparisons = function(data, ..., xlim=NULL, ylim=NULL, title = "", legend.position=c(0.9,0.1), legend.key.width = unit(25, "mm"), reference=NULL) { + # If no reference mortality table is given, use the first table (data if its a valuation table) + if (missing(reference)) { + if (inherits(data, "valuationTable")) { + reference = data; + } else { + reference = NULL;# TODO; + } + } + refname=reference@name; + yAxisLabel = expression(paste("Sterbewahrscheinlichkeit ", q[x], " relativ zu ", refname)); +yAxisLabel[1][3]="asfd" +yAxisLabel + if (!is.data.frame(data)) { + data = makeQxDataFrame(data, ..., reference=reference); + } + + pl = ggplot(data, aes(x = x, y = y, colour = data$group)) + + theme_bw() + + theme( + plot.title = element_text(size=18, face="bold"), + legend.title = element_text(size=14, face="bold.italic"), + # legend in bottom right corner of the plot + legend.justification=c(1,0), legend.position=legend.position, + # No box around legend entries + legend.key = element_blank(), + legend.key.width = legend.key.width, + legend.background = element_rect(colour="gray50", linetype="solid") + ) + + geom_line() + + coord_cartesian(xlim=xlim, ylim=ylim) + + scale_y_continuous( +# #name=substitute(paste("Sterbewahrscheinlichkeit ", q[x], " relativ zu ", refname), env=list(refname=refname))#, + name=substitute(refname^3, env=list(refname=refname)), + labels=percent +# # breaks = scales::trans_breaks('log10', function(x) 10^x), +# # labels = scales::trans_format('log10', scales::math_format(10^.x)) +# #minor_breaks = log(c(sapply(x, function(x) seq(0, x, x/10))), 10) + ) + + scale_x_continuous( + name="Alter", + #breaks = function (limits) scales::trans_breaks('', function(x) 10^x), + breaks = function (limits) seq(max(min(limits),0),max(limits),5), + minor_breaks = function (limits) seq(max(round(min(limits)),0),round(max(limits)),1), + #labels = scales::trans_format('log10', scales::math_format(10^.x)) + + ) + + # annotation_logticks(sides="lr") + + xlab("Alter") + labs(colour="Sterbetafel"); + if (title != "") { + pl = pl + ggtitle(title); + } + pl +} + # # plotValuationTables(mort.AT.census.1869.male, mort.AT.census.1869.female, mort.AT.census.2011.male, mort.AT.census.2011.female, AVOe2005R.male, AVOe2005R.female, YOB=1972,title="Vergleich österreichische Sterbetafeln, YOB=1972 (bei Generationentafeln)") #