Commit 28a65938 authored by Kevin's avatar Kevin

added median bias

parent 101dd90d
...@@ -77,7 +77,7 @@ rule generate_validation_data: ...@@ -77,7 +77,7 @@ rule generate_validation_data:
# adjust threads by model type # adjust threads by model type
def get_rule_threads(wildcards): def get_rule_threads(wildcards):
if wildcards.model in ("locf", "msm"): if wildcards.model in ("locf", "msm", "msm_age"):
return 1 return 1
else: else:
return config["stan"]["chains"] return config["stan"]["chains"]
...@@ -98,17 +98,6 @@ rule fit_model_validation_set: ...@@ -98,17 +98,6 @@ rule fit_model_validation_set:
# helper rule to just build all posterior datasets
rule model_posteriors:
input:
["output/v1.1/data/validation/posteriors/%s/df_posterior_mi_%i_fold_%i.rds" % (m, i, j)
for m in ("locf", "msm", "msm_age", "gp", "gp_nb", "mm", "mm_nb")
for i in range(1, config["mi_m"] + 1)
for j in range(1, config["folds"] + 1)
]
# rules for imputing on entire dataset # rules for imputing on entire dataset
rule generate_imputation_data: rule generate_imputation_data:
input: input:
...@@ -181,22 +170,26 @@ rule impute_msm_v1_1: ...@@ -181,22 +170,26 @@ rule impute_msm_v1_1:
# helper rule to just build all posterior datasets
rule model_posteriors_v1_1:
input:
["output/v1.1/data/validation/posteriors/%s/df_posterior_mi_%i_fold_%i.rds" % (m, i, j)
for m in ("locf", "msm", "msm_age", "gp", "gp_nb", "mm", "mm_nb")
for i in range(1, config["mi_m"] + 1)
for j in range(1, config["folds"] + 1)
]
rule manuscript: rule manuscript_v1_1:
input: input:
pop_report = rules.prepare_data.output, pop_report = rules.prepare_data.output,
posteriors = rules.model_posteriors.input, posteriors = rules.model_posteriors.input,
markdown = "manuscript/manuscript.Rmd" markdown = "manuscript/manuscript.Rmd"
output: output:
pdf = "output/{version}/manuscript.docx", pdf = "output/v1.1/manuscript.docx",
figures = "output/{version}/manuscript_figures.zip" figures = "output/v1.1/manuscript_figures.zip"
shell: shell:
""" """
mkdir -p output/{wildcards.version} mkdir -p output/{wildcards.version}
Rscript -e "rmarkdown::render(\\"{input.markdown}\\", output_dir = \\"output/{wildcards.version}\\", params = list(data_dir = \\"../output/{wildcards.version}/data\\", config_file = \\"../config.yml\\"))" Rscript -e "rmarkdown::render(\\"{input.markdown}\\", output_dir = \\"output/v1.1\\", params = list(data_dir = \\"../output/v1.1/data\\", config_file = \\"../config.yml\\"))"
mv manuscript/figures.zip {output.figures} mv manuscript/figures.zip {output.figures}
""" """
rule manuscript_v1_1:
input:
"output/v1.1/manuscript.docx"
...@@ -545,14 +545,15 @@ p1 <- df_predictions %>% ...@@ -545,14 +545,15 @@ p1 <- df_predictions %>%
summarize( summarize(
RMSE = mean((GOSE - prediction)^2, na.rm = TRUE) %>% sqrt, RMSE = mean((GOSE - prediction)^2, na.rm = TRUE) %>% sqrt,
MAE = mean(abs(GOSE - prediction), na.rm = TRUE), MAE = mean(abs(GOSE - prediction), na.rm = TRUE),
bias = mean(prediction, na.rm = TRUE) - mean(GOSE, na.rm = TRUE) Bias = mean(prediction, na.rm = TRUE) - mean(GOSE, na.rm = TRUE),
`Pr[est > true] - Pr[est < true]` = mean(prediction > GOSE, na.rm = TRUE) - mean(prediction < GOSE, na.rm = TRUE)
) %>% ) %>%
ungroup %>% ungroup %>%
gather(error, value, -model, -fold) %>% gather(error, value, -model, -fold) %>%
ggplot(aes(model, value)) + ggplot(aes(model, value)) +
geom_hline(yintercept = 0, color = "black") + geom_hline(yintercept = 0, color = "black") +
geom_boxplot() + geom_boxplot() +
facet_wrap(~error) + facet_wrap(~error, nrow = 1) +
scale_y_continuous(name = "", breaks = seq(-2, 8, .25), limits = c(-.5, 1.5)) + scale_y_continuous(name = "", breaks = seq(-2, 8, .25), limits = c(-.5, 1.5)) +
scale_x_discrete("") + scale_x_discrete("") +
theme_bw() + theme_bw() +
...@@ -569,14 +570,15 @@ p2 <- df_predictions %>% ...@@ -569,14 +570,15 @@ p2 <- df_predictions %>%
summarize( summarize(
RMSE = mean((GOSE - prediction)^2, na.rm = TRUE) %>% sqrt, RMSE = mean((GOSE - prediction)^2, na.rm = TRUE) %>% sqrt,
MAE = mean(abs(GOSE - prediction), na.rm = TRUE), MAE = mean(abs(GOSE - prediction), na.rm = TRUE),
bias = mean(prediction, na.rm = TRUE) - mean(GOSE, na.rm = TRUE) Bias = mean(prediction, na.rm = TRUE) - mean(GOSE, na.rm = TRUE),
`Pr[est > true] - Pr[est < true]` = mean(prediction > GOSE, na.rm = TRUE) - mean(prediction < GOSE, na.rm = TRUE)
) %>% ) %>%
ungroup %>% ungroup %>%
gather(error, value, -model, -fold) %>% gather(error, value, -model, -fold) %>%
ggplot(aes(model, value)) + ggplot(aes(model, value)) +
geom_hline(yintercept = 0, color = "black") + geom_hline(yintercept = 0, color = "black") +
geom_boxplot() + geom_boxplot() +
facet_wrap(~error) + facet_wrap(~error, nrow = 1) +
scale_y_continuous(name = "", breaks = seq(-2, 8, .25), limits = c(-.5, 1.5)) + scale_y_continuous(name = "", breaks = seq(-2, 8, .25), limits = c(-.5, 1.5)) +
scale_x_discrete("") + scale_x_discrete("") +
theme_bw() + theme_bw() +
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment