Revision 2d03921e
Added by Benoit Parmentier about 11 years ago
| climate/research/oregon/interpolation/analyses_papers_methods_comparison_part5.R | ||
|---|---|---|
| 59 | 59 |
return(prop_obj) |
| 60 | 60 |
} |
| 61 | 61 |
|
| 62 |
#Calculate the difference between training and testing in two different data.frames. Columns to substract are provided. |
|
| 63 |
diff_df<-function(tb_s,tb_v,list_metric_names){
|
|
| 64 |
tb_diff<-vector("list", length(list_metric_names))
|
|
| 65 |
for (i in 1:length(list_metric_names)){
|
|
| 66 |
metric_name<-list_metric_names[i] |
|
| 67 |
tb_diff[[i]] <-tb_s[,c(metric_name)] - tb_v[,c(metric_name)] |
|
| 68 |
} |
|
| 69 |
names(tb_diff)<-list_metric_names |
|
| 70 |
tb_diff<-as.data.frame(do.call(cbind,tb_diff)) |
|
| 71 |
return(tb_diff) |
|
| 72 |
} |
|
| 73 |
|
|
| 62 | 74 |
################## PARAMETERS ########## |
| 63 | 75 |
|
| 64 | 76 |
|
| ... | ... | |
| 80 | 92 |
in_dir13 <-"/data/project/layers/commons/Oregon_interpolation/output_data_365d_gwr_cai_lst_comb3_09282013" |
| 81 | 93 |
in_dir14 <-"/data/project/layers/commons/Oregon_interpolation/output_data_365d_gwr_fus_lst_comb3_09232013" |
| 82 | 94 |
in_dir15 <-"/data/project/layers/commons/Oregon_interpolation/output_data_365d_gwr_fus_lst_comb3_09262013" |
| 95 |
in_dir16 <-"/data/project/layers/commons/Oregon_interpolation/output_data_365d_gam_cai_lst_comb3_10042013" |
|
| 96 |
|
|
| 83 | 97 |
|
| 84 | 98 |
#better as list and load one by one specific element from the object |
| 85 | 99 |
raster_prediction_obj1 <-load_obj(file.path(in_dir1,"raster_prediction_obj_gam_CAI_dailyTmax_365d_gam_CAI_lst_comb3_08312013.RData")) |
| ... | ... | |
| 99 | 113 |
raster_prediction_obj14 <-load_obj(file.path(in_dir14,"raster_prediction_obj_gwr_fusion_dailyTmax_365d_gwr_fus_lst_comb3_09232013.RData")) |
| 100 | 114 |
raster_prediction_obj15 <-load_obj(file.path(in_dir15,"raster_prediction_obj_gwr_fusion_dailyTmax_365d_gwr_fus_lst_comb3_09262013.RData")) |
| 101 | 115 |
|
| 102 |
#raster_prediction_obj_gwr_CAI_dailyTmax_365d_gwr_cai_lst_comb3_09282013.RData |
|
| 103 |
#raster_prediction_obj_gwr_fusion_dailyTmax_365d_gwr_fus_lst_comb3_09262013.RData |
|
| 116 |
raster_prediction_obj16 <-load_obj(file.path(in_dir16,"raster_prediction_obj_gam_CAI_dailyTmax_365d_gam_cai_lst_comb3_10042013.RData")) |
|
| 104 | 117 |
|
| 105 | 118 |
out_dir<-"/home/parmentier/Data/IPLANT_project/paper_multitime_scale__analyses_tables_fig_09032013" |
| 106 | 119 |
setwd(out_dir) |
| 107 | 120 |
y_var_name <- "dailyTmax" |
| 108 | 121 |
y_var_month <- "TMax" |
| 109 | 122 |
#y_var_month <- "LSTD_bias" |
| 110 |
out_suffix <- "_OR_09292013"
|
|
| 123 |
out_suffix <- "_OR_10102013"
|
|
| 111 | 124 |
#script_path<-"/data/project/layers/commons/data_workflow/env_layers_scripts/" |
| 112 | 125 |
#### FUNCTION USED IN SCRIPT |
| 113 | 126 |
|
| ... | ... | |
| 217 | 230 |
|
| 218 | 231 |
##### Calculate differences |
| 219 | 232 |
|
| 220 |
#Calculate the difference between training and testing in two different data.frames. Columns to substract are provided. |
|
| 221 |
diff_df<-function(tb_s,tb_v,list_metric_names){
|
|
| 222 |
tb_diff<-vector("list", length(list_metric_names))
|
|
| 223 |
for (i in 1:length(list_metric_names)){
|
|
| 224 |
metric_name<-list_metric_names[i] |
|
| 225 |
tb_diff[[i]] <-tb_s[,c(metric_name)] - tb_v[,c(metric_name)] |
|
| 226 |
} |
|
| 227 |
names(tb_diff)<-list_metric_names |
|
| 228 |
tb_diff<-as.data.frame(do.call(cbind,tb_diff)) |
|
| 229 |
return(tb_diff) |
|
| 230 |
} |
|
| 231 |
|
|
| 232 |
|
|
| 233 | 233 |
metric_names <- c("mae","rmse","me","r")
|
| 234 | 234 |
diff_kriging_CAI <- diff_df(tb_s_kriging_CAI,tb_v_kriging_CAI,metric_names) |
| 235 |
diff_gam_CAI <- diff_df(tb_s_gam_CAI,tb_v_gam_CAI,metric_names) |
|
| 235 |
diff_gam_CAI <- diff_df(tb_s_gam_CAI[tb_s_gam_CAI$pred_mod!="mod_kr"],tb_v_gam_CAI,metric_names)
|
|
| 236 | 236 |
diff_gwr_CAI <- diff_df(tb_s_gwr_CAI,tb_v_gwr_CAI,metric_names) |
| 237 | 237 |
|
| 238 | 238 |
layout_m<-c(1,1) #one row two columns |
| ... | ... | |
| 301 | 301 |
|
| 302 | 302 |
### NOW PLOT OF COMPARISON BETWEEN Kriging and GAM |
| 303 | 303 |
|
| 304 |
#Now get |
|
| 304 |
#Now get variance and range for holdout an dmethods.
|
|
| 305 | 305 |
|
| 306 | 306 |
tb_v_gam_CAI |
| 307 | 307 |
tb_v_gam_fus |
| ... | ... | |
| 342 | 342 |
test2<-test[test$method_interp%in% c("gam_fus","gam_CAI"),]
|
| 343 | 343 |
test2[1:24,] |
| 344 | 344 |
#head(ac_prop_tb) |
| 345 |
test3<-subset(test,prop_month==0 & method_interp%in%c("gam_CAI"))
|
|
| 346 |
#test3<-test[test$method_interp%in% c("gam_CAI"),]
|
|
| 347 |
test3[1:24,] |
|
| 345 | 348 |
|
| 346 | 349 |
#list_prop_obj$avg_tb |
| 347 | 350 |
#xyplot(as.formula(plot_formula),group=pred_mod,type="b", |
| ... | ... | |
| 352 | 355 |
# pch=1:length(avg_tb$pred_mod))), |
| 353 | 356 |
# auto.key=list(columns=5)) |
| 354 | 357 |
|
| 358 |
## DAILY DEVIATIONS WITH MODELS |
|
| 359 |
### Examining results with models for daily deviation called 2: need to subset the models |
|
| 360 |
|
|
| 361 |
#c("mod1.dev_mod1","mod1.dev_mod2","mod2.dev_mod1","mod2.dev_mod2","mod3.dev_mod1","mod3.dev_mod2","mod4.dev_mod1",
|
|
| 362 |
# "mod4.dev_mod2","mod5.dev_mod1","mod5.dev_mod2","mod6.dev_mod1","mod6.dev_mod2","mod7.dev_mod1","mod7.dev_mod2", |
|
| 363 |
# "mod_kr.dev_mod1" "mod_kr.dev_mod2") |
|
| 364 |
|
|
| 365 |
tb_s_gam_CAI_selected <-subset(tb_s_gam_CAI,prop_month==0 & pred_mod%in%c("mod1","mod4","mod5","mod6","mod7","mod8"))
|
|
| 366 |
tb_v_gam_CAI_selected <-subset(tb_v_gam_CAI,prop_month==0 & pred_mod%in%c("mod1","mod4","mod5","mod6","mod7","mod8"))
|
|
| 367 |
tb_s_gwr_CAI_selected <-subset(tb_s_gwr_CAI,prop_month==0 & pred_mod%in%c("mod1","mod4","mod5","mod6","mod7","mod8"))
|
|
| 368 |
tb_v_gwr_CAI_selected <-subset(tb_v_gwr_CAI,prop_month==0 & pred_mod%in%c("mod1","mod4","mod5","mod6","mod7","mod8"))
|
|
| 369 |
tb_s_kriging_CAI_selected <-subset(tb_s_kriging_CAI,prop_month==0 & pred_mod%in%c("mod1","mod4","mod5","mod6","mod7","mod8"))
|
|
| 370 |
tb_v_kriging_CAI_selected <-subset(tb_v_kriging_CAI,prop_month==0 & pred_mod%in%c("mod1","mod4","mod5","mod6","mod7","mod8"))
|
|
| 371 |
|
|
| 372 |
tb_s_gam_CAI2_selected <- subset(raster_prediction_obj16$tb_diagnostic_s,pred_mod%in%c("mod1.dev_mod1","mod1.dev_mod2","mod3.dev_mod1",
|
|
| 373 |
"mod3.dev_mod2","mod4.dev_mod1", "mod4.dev_mod2", |
|
| 374 |
"mod5.dev_mod1","mod5.dev_mod2","mod6.dev_mod1", |
|
| 375 |
"mod6.dev_mod2","mod7.dev_mod1","mod7.dev_mod2")) |
|
| 376 |
tb_v_gam_CAI2_selected <- subset(raster_prediction_obj16$tb_diagnostic_v,pred_mod%in%c("mod1.dev_mod1","mod1.dev_mod2","mod3.dev_mod1",
|
|
| 377 |
"mod3.dev_mod2","mod4.dev_mod1", "mod4.dev_mod2", |
|
| 378 |
"mod5.dev_mod1","mod5.dev_mod2","mod6.dev_mod1", |
|
| 379 |
"mod6.dev_mod2","mod7.dev_mod1","mod7.dev_mod2")) |
|
| 380 |
|
|
| 381 |
metric_names <- c("mae","rmse","me","r")
|
|
| 382 |
|
|
| 383 |
diff_gam_cai2 <- diff_df(tb_s_gam_CAI2_selected,tb_v_gam_CAI2_selected,metric_names) |
|
| 384 |
diff_kriging_CAI_selected <- diff_df(tb_s_kriging_CAI_selected,tb_v_kriging_CAI_selected,metric_names) |
|
| 385 |
diff_gam_CAI_selected <- diff_df(tb_s_gam_CAI_selected,tb_v_gam_CAI_selected,metric_names) |
|
| 386 |
diff_gwr_CAI_selected <- diff_df(tb_s_gwr_CAI_selected,tb_v_gwr_CAI_selected,metric_names) |
|
| 387 |
|
|
| 388 |
layout_m<-c(1,1) #one row two columns |
|
| 389 |
par(mfrow=layout_m) |
|
| 390 |
|
|
| 391 |
png(paste("Figure__accuracy_rmse_prop_month_",plot_name,out_suffix,".png", sep=""),
|
|
| 392 |
height=480*layout_m[1],width=480*layout_m[2]) |
|
| 393 |
boxplot(diff_gam_cai2$rmse,diff_gam_CAI_selected$rmse, |
|
| 394 |
diff_kriging_CAI_selected$rmse,diff_gwr_CAI_selected$rmse,names=c("gam_cai2","gam_CAI","kriging","gwr"))
|
|
| 395 |
|
|
| 396 |
dev.off() |
|
Also available in: Unified diff
analyses multi-time scale part 5, examining results with covariates in daily deviations