From da498de171f06bd469e149324904afe2d46e09b8 Mon Sep 17 00:00:00 2001
From: qianyu313 <qianyu313@gmail.com>
Date: Wed, 14 Feb 2024 11:30:54 -0800
Subject: [PATCH] interval plot for national and admin 1 results

---
 R/intervalPlot.R    | 222 ++++++++++++++++++++++++++++++++++++++------
 man/intervalPlot.Rd |  12 ++-
 2 files changed, 202 insertions(+), 32 deletions(-)

diff --git a/R/intervalPlot.R b/R/intervalPlot.R
index f6c3dc0..5157f55 100644
--- a/R/intervalPlot.R
+++ b/R/intervalPlot.R
@@ -2,7 +2,9 @@
 #'
 #' This function return scatter plot at admin 1 level for any two model results
 #'
-#' @param res  model result  using surveyPrev
+#' @param model  list of model results  using surveyPrev
+#' @param admin  level of plot
+#' @param  group plot by group or not
 #' @return This function returns the dataset that contain district name and population for given  tiff files and polygons of admin level.
 #' @import ggplot2
 #' @author Qianyu Dong
@@ -61,44 +63,204 @@
 #' }
 #'
 #' @export
+#'
+
+
+intervalPlot <- function(admin= 0, model=list( list("name 1"= fit1, "name 1"= fit2)), group=F){
+
+
+
+
+  if(admin==0){
+
+    dt=data.frame(mean=NA, lower=NA, upper=NA, model=NA,group=NA)
+
+
+
+
+    for (i in 1:length(model)) {
+
+      if(colnames(model[[i]]$agg.natl[1])=="direct.est"){
+        colnames(model[[i]]$agg.natl)[colnames(model[[i]]$agg.natl) == 'direct.est'] <- 'mean'
+        colnames(model[[i]]$agg.natl)[colnames(model[[i]]$agg.natl) == 'direct.lower'] <- 'lower'
+        colnames(model[[i]]$agg.natl)[colnames(model[[i]]$agg.natl) == 'direct.upper'] <- 'upper'
+      }
+
+      model[[i]]$agg.natl$model= names(model[i])
+      dt[i,]=model[[i]]$agg.natl[,c("mean","lower","upper","model")]
+    }
+
+
+    if(group){
+      for (i in 1:length(model)) {
+        dt[i,]$group= model[[i]]$group
+      }
+    }
 
-intervalPlot <- function(res){
+    rownames(dt)<-dt$model
 
-  # TODO: refer by variable name and make it work with both admin 1 and 2 cases
-  data<-res[[1]]
-  linedata<-res[[2]]
-  line2<-res[[3]]$mean
 
-  plot_fun <- function(dat) {
-    line1 = linedata[unique(dat$admin1.name),"mean"]
-    g <- ggplot(dat)
+    if(group){
+      ggplot(dt, aes(y = model, x = mean,group = group)) +
+        geom_point(aes(shape=group)) +
+        geom_errorbarh(aes(xmin = lower , xmax = upper), alpha = 0.5)
 
-    if("type" %in% colnames(dat)){
-      g <- g + aes(x = admin2.name, y = mean, color = type)
     }else{
-      g <- g + aes(x = admin2.name, y = mean)
+      ggplot(dt, aes(y = model, x = mean)) +
+        geom_point(aes()) +
+        geom_errorbarh(aes(xmin = lower , xmax = upper), alpha = 0.5)
+
     }
-    g <- g +
-        geom_point( position = position_dodge(width = 0.5),size = .8) +
-        geom_hline(  aes(yintercept =line1,linetype = "dashed"),color="#d95f02",linewidth = .8) +
-        geom_hline( aes(yintercept =line2, linetype = "solid"),color="#d95f02",linewidth = .8) +
-        geom_errorbar(aes(ymin = lower , ymax = upper), alpha = 0.7,position = position_dodge(width = 0.5), width = 0.2) +
-        scale_color_brewer(palette = "Set1") +
-      labs(title = unique(dat$admin1.name)) +
-      xlab("") + ylab("") +
-      scale_linetype_manual(values = c("dashed", "solid"),
-                            labels = c("admin1","national")) +
-      theme_bw() +
-      theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))
+
+
+
+
+  }else if (admin==1){
+
+    dt=data.frame(admin1.name=rep(NA,500),mean=rep(NA,500), lower=rep(NA,500), upper=rep(NA,500), model=rep(NA,500),group=rep(NA,500))
+
+    for (i in 1:length(model)) {
+
+
+
+
+      if( is.null(model[[i]]$agg.admin1)){
+
+        if(colnames(model[[i]]$res.admin1)[2]=="direct.est"){
+          colnames(model[[i]]$res.admin1)[colnames(model[[i]]$res.admin1) == 'direct.est'] <- 'mean'
+          colnames(model[[i]]$res.admin1)[colnames(model[[i]]$res.admin1) == 'direct.lower'] <- 'lower'
+          colnames(model[[i]]$res.admin1)[colnames(model[[i]]$res.admin1) == 'direct.upper'] <- 'upper'
+        }
+        model[[i]]$res.admin1$model= names(model[i])
+
+        if( !is.null(model[[i]]$res.admin1$type)){
+          model[[i]]$res.admin1=model[[i]]$res.admin1[model[[i]]$res.admin1$type=="full",]
+        }
+
+
+        dd=dim(model[[i]]$res.admin1)[1]
+
+
+        dt[(1+ (i-1)*dd):(i*dd),]= model[[i]]$res.admin1[,c("admin1.name","mean","lower","upper","model")]
+        if(group){
+          dt[(1+ (i-1)*dd):(i*dd),]$group= model[[i]]$group
+        }
+      }else{
+
+        if(colnames(model[[i]]$agg.admin1)[2]=="direct.est"){
+          colnames(model[[i]]$agg.admin1)[colnames(model[[i]]$agg.admin1) == 'direct.est'] <- 'mean'
+          colnames(model[[i]]$agg.admin1)[colnames(model[[i]]$agg.admin1) == 'direct.lower'] <- 'lower'
+          colnames(model[[i]]$agg.admin1)[colnames(model[[i]]$agg.admin1) == 'direct.upper'] <- 'upper'
+        }
+        model[[i]]$agg.admin1$model= names(model[i])
+
+        if( !is.null(model[[i]]$agg.admin1$type)){
+          model[[i]]$agg.admin1=model[[i]]$agg.admin1[model[[i]]$agg.admin1$type=="full",]
+        }
+
+
+        dd=dim(model[[i]]$agg.admin1)[1]
+        dt[(1+ (i-1)*dd):(i*dd),]= model[[i]]$agg.admin1[,c("admin1.name","mean","lower","upper","model")]
+        if(group){
+          dt[(1+ (i-1)*dd):(i*dd),]$group= model[[i]]$group
+        }
+      }
+    }
+
+    dt=na.omit(dt)
+
+
+
+
+    if(group){
+      ggplot(dt, aes(x = admin1.name, y = mean, group = model, color = model,shape=group)) +
+        geom_point( position = position_dodge(width = 0.7)) +
+        scale_shape_manual(values = c(0:5, 15:19)) +
+        geom_errorbar(aes(ymin = lower,
+                          ymax = upper, group = model), alpha = 0.8, position = position_dodge(width = 0.7))+
+        # scale_color_brewer(palette="Set3") +
+        theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
+        labs(title = "", x = "Region", y = "value")+
+        theme(legend.title = element_text(size=10), # Increase legend title size
+              legend.text = element_text(size=12), # Increase legend text size
+              legend.key.size = unit(1.5, 'lines'), # Increase legend key size
+              axis.text.x = element_text(size=12), # Increase x axis text size
+              axis.text.y = element_text(size=12))
+
+    }else{
+      ggplot(dt, aes(x = admin1.name, y = mean, group = model, color = model)) +
+        geom_point( position = position_dodge(width = 0.7)) +
+        scale_shape_manual(values = c(0:5, 15:19)) +
+        geom_errorbar(aes(ymin = lower,
+                          ymax = upper, group = model), alpha = 0.8, position = position_dodge(width = 0.7))+
+        # scale_color_brewer(palette="Set3") +
+        theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
+        labs(title = "", x = "Region", y = "value")+
+        theme(legend.title = element_text(size=10), # Increase legend title size
+              legend.text = element_text(size=12), # Increase legend text size
+              legend.key.size = unit(1.5, 'lines'), # Increase legend key size
+              axis.text.x = element_text(size=12), # Increase x axis text size
+              axis.text.y = element_text(size=12))
+
+    }
+
+  }else{
+
+
+
+
+
+
+
   }
 
- plots <- NULL
- for (adm1 in unique(data$admin1.name)){
-  g <- plot_fun(subset(data, admin1.name == adm1))
-  plots[[adm1]] <- g
- }
 
-  return(plots)
+
+
+
+
+
+
+
+
+
+
+ #
+ #
+ #  data<-res[[1]]
+ #  linedata<-res[[2]]
+ #  line2<-res[[3]]$mean
+ #
+ #  plot_fun <- function(dat) {
+ #    line1 = linedata[unique(dat$admin1.name),"mean"]
+ #    g <- ggplot(dat)
+ #
+ #    if("type" %in% colnames(dat)){
+ #      g <- g + aes(x = admin2.name, y = mean, color = type)
+ #    }else{
+ #      g <- g + aes(x = admin2.name, y = mean)
+ #    }
+ #    g <- g +
+ #        geom_point( position = position_dodge(width = 0.5),size = .8) +
+ #        geom_hline(  aes(yintercept =line1,linetype = "dashed"),color="#d95f02",linewidth = .8) +
+ #        geom_hline( aes(yintercept =line2, linetype = "solid"),color="#d95f02",linewidth = .8) +
+ #        geom_errorbar(aes(ymin = lower , ymax = upper), alpha = 0.7,position = position_dodge(width = 0.5), width = 0.2) +
+ #        scale_color_brewer(palette = "Set1") +
+ #      labs(title = unique(dat$admin1.name)) +
+ #      xlab("") + ylab("") +
+ #      scale_linetype_manual(values = c("dashed", "solid"),
+ #                            labels = c("admin1","national")) +
+ #      theme_bw() +
+ #      theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))
+ #  }
+ #
+ # plots <- NULL
+ # for (adm1 in unique(data$admin1.name)){
+ #  g <- plot_fun(subset(data, admin1.name == adm1))
+ #  plots[[adm1]] <- g
+ # }
+ #
+ #  return(plots)
 
 }
 
diff --git a/man/intervalPlot.Rd b/man/intervalPlot.Rd
index b4f169b..389df02 100644
--- a/man/intervalPlot.Rd
+++ b/man/intervalPlot.Rd
@@ -4,10 +4,18 @@
 \alias{intervalPlot}
 \title{Get scatter plot for any two model results}
 \usage{
-intervalPlot(res)
+intervalPlot(
+  admin = 0,
+  model = list(list(`name 1` = fit1, `name 1` = fit2)),
+  group = F
+)
 }
 \arguments{
-\item{res}{model result  using surveyPrev}
+\item{admin}{level of plot}
+
+\item{model}{list of model results  using surveyPrev}
+
+\item{group}{plot by group or not}
 }
 \value{
 This function returns the dataset that contain district name and population for given  tiff files and polygons of admin level.