---
title: "02_Visualization"
output: html_document
date: "2023-11-07"
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
```{r Load data and packages}
library(ggplot2);library(ggpubr)
regression_table <- read.csv("Outputs/regression_table.csv")
regression_table_total <- read.csv("Outputs/regression_table_total.csv")
regression_table_total_commune <- read.csv("Outputs/regression_table_total_commune.csv")
regression_table_3km <- read.csv("Outputs/regression_table_3km.csv")
regression_table_6km <- read.csv("Outputs/regression_table_6km.csv")
regression_table_9km <- read.csv("Outputs/regression_table_9km.csv")
sensitivity_table <- read.csv("Outputs/sensitivity_table.csv")
agg_3km <- read.csv("Outputs/regression_table_3km_agg.csv")
agg_combined <- read.csv("Outputs/regression_table_combined_agg.csv")
```
```{r Main analysis plot}
regression_table_main <- rbind(regression_table, regression_table_3km, regression_table_total, regression_table_total_commune)
regression_table_main$Lower <- regression_table_main$Estimate - 1.96 * regression_table_main$Cluster.s.e.
regression_table_main$Upper <- regression_table_main$Estimate + 1.96 * regression_table_main$Cluster.s.e.
regression_table_main$Lower90 <- regression_table_main$Estimate - qnorm(0.95) * regression_table_main$Cluster.s.e.
regression_table_main$Upper90 <- regression_table_main$Estimate + qnorm(0.95) * regression_table_main$Cluster.s.e.
regression_table_main$Analysis[regression_table_main$Analysis=="Commune"]<- "Entrance"
regression_table_main$Analysis[regression_table_main$Analysis=="PA"] <- "Total PA"
regression_table_main$Location[regression_table_main$Location=="Inside PA + Buffer"]<- "Inside PA + Buffer 3km"
regression_table_main$Analysis <- factor(regression_table_main$Analysis, levels=c("Total PA", "Entrance"))
regression_table_main$Location <- factor(regression_table_main$Location, levels=c("Inside PA", "Buffer 3km", "Inside PA + Buffer 3km"))
main_plot <- ggplot(data=regression_table_main[regression_table_main$Variable=="Tourism",], aes(x = Location , y= Estimate, color=Analysis))+
geom_hline(aes(yintercept=0), linetype="dashed")+
geom_linerange(aes(ymin = Lower, ymax = Upper), position = position_dodge(0.3))+
geom_linerange(aes(ymin = Lower90, ymax = Upper90), position = position_dodge(0.3), linewidth=1.5)+
geom_point(aes(x = Location , y= Estimate), position = position_dodge(0.3), size=3.5)+
theme_classic()+
scale_color_manual(values=c("magenta4", "thistle3"))+
xlab("")+
ylab("Estimated Impact of Tourism on Deforestation")
write.csv(regression_table_main, "Outputs/regression_table_MAIN.csv")
```
```{r Sensitivity analysis/ robustness check plot}
sensitivity_table$Location <- "Inside PA"
sensitivity_table$Lower <- sensitivity_table$Estimate - 1.96 * sensitivity_table$Cluster.s.e.
sensitivity_table$Upper <- sensitivity_table$Estimate + 1.96 * sensitivity_table$Cluster.s.e.
sensitivity_table$Lower90 <- sensitivity_table$Estimate - qnorm(0.95) * sensitivity_table$Cluster.s.e.
sensitivity_table$Upper90 <- sensitivity_table$Estimate + qnorm(0.95) * sensitivity_table$Cluster.s.e.
#Distance from PA boundary
regression_table_distance <- rbind( regression_table_6km,regression_table_9km )
regression_table_distance$Lower <- regression_table_distance$Estimate - 1.96 * regression_table_distance$Cluster.s.e.
regression_table_distance$Upper <- regression_table_distance$Estimate + 1.96 * regression_table_distance$Cluster.s.e.
regression_table_distance$Lower90 <- regression_table_distance$Estimate - qnorm(0.95) * regression_table_distance$Cluster.s.e.
regression_table_distance$Upper90 <- regression_table_distance$Estimate + qnorm(0.95)* regression_table_distance$Cluster.s.e.
regression_table_distance$Location <- factor(regression_table_distance$Location, levels=c("Buffer 6km", "Buffer 9km"))
regression_table_distance$Analysis[regression_table_distance$Analysis=="Commune"]<- "Entrance"
regression_table_distance$Analysis[regression_table_distance$Analysis=="PA"]<- "Total PA"
regression_table_distance$Analysis <- factor(regression_table_distance$Analysis, levels=c("Total PA", "Entrance"))
distance_plot <- ggplot(data=regression_table_distance[regression_table_distance$Variable=="Tourism",], aes(x = Location , y= Estimate, color=Analysis))+
geom_hline(aes(yintercept=0), linetype="dashed")+
geom_linerange(aes(ymin = Lower, ymax = Upper), position = position_dodge(0.3))+
geom_linerange(aes(ymin = Lower90, ymax = Upper90), position = position_dodge(0.3), linewidth=1.5)+
geom_point(position = position_dodge(0.3), size=3.5)+
theme_classic()+
scale_color_manual(values=c("magenta4", "thistle3"))+
theme(legend.position="top")+
ylab("Estimated Impact of Tourism on Deforestation")
write.csv(regression_table_distance, "Outputs/regression_table_DISTANCE.csv")
#Aggregate Analysis
sensitivity_table_aggregate <- sensitivity_table[sensitivity_table$model=="Aggregate",]
agg_3km$Lower <- agg_3km$Estimate - 1.96 * agg_3km$Cluster.s.e.
agg_3km$Upper <- agg_3km$Estimate + 1.96 * agg_3km$Cluster.s.e.
agg_combined$Lower <- agg_combined$Estimate - 1.96 * agg_combined$Cluster.s.e.
agg_combined$Upper <- agg_combined$Estimate + 1.96 * agg_combined$Cluster.s.e.
agg_3km$Lower90 <- agg_3km$Estimate - qnorm(0.95) * agg_3km$Cluster.s.e.
agg_3km$Upper90 <- agg_3km$Estimate + qnorm(0.95) * agg_3km$Cluster.s.e.
agg_combined$Lower90 <- agg_combined$Estimate - qnorm(0.95) * agg_combined$Cluster.s.e.
agg_combined$Upper90 <- agg_combined$Estimate + qnorm(0.95)* agg_combined$Cluster.s.e.
colnames(sensitivity_table_aggregate)[6]<- "Analysis"
agg_total <- rbind(sensitivity_table_aggregate, agg_3km, agg_combined)
agg_total$Location[agg_total$Location=="Inside PA + Buffer"] <- "Inside PA + Buffer 3km"
agg_total$Location <- factor(agg_total$Location, levels= c("Inside PA", "Buffer 3km", "Inside PA + Buffer 3km"))
aggregate_plot <- ggplot(data=agg_total[agg_total$Variable=="Tourism",], aes(x = Location , y= Estimate, color=Analysis))+
geom_hline(aes(yintercept=0), linetype="dashed")+
geom_linerange(aes(ymin = Lower, ymax = Upper), position = position_dodge(0.3), color="magenta4")+
geom_linerange(aes(ymin = Lower90, ymax = Upper90), position = position_dodge(0.3), color="magenta4", linewidth=1.5)+
geom_point(position = position_dodge(0.3), size=3.5, color="magenta4")+
theme_classic()+
theme(legend.position="top")+
ylab("Estimated Impact of Tourism on Deforestation")
write.csv(agg_total, "Outputs/regression_table_AGGREGATE.csv")
#Rival explanations
sensitivity_table_rival_buffer <- read.csv("Outputs/regression_table_3km_precip.csv")
sensitivity_table_rival_buffer$Lower <- sensitivity_table_rival_buffer$Estimate - 1.96 * sensitivity_table_rival_buffer$Cluster.s.e.
sensitivity_table_rival_buffer$Upper <- sensitivity_table_rival_buffer$Estimate + 1.96 * sensitivity_table_rival_buffer$Cluster.s.e.
sensitivity_table_rival_buffer$Lower90 <- sensitivity_table_rival_buffer$Estimate - qnorm(0.95) * sensitivity_table_rival_buffer$Cluster.s.e.
sensitivity_table_rival_buffer$Upper90 <- sensitivity_table_rival_buffer$Estimate + qnorm(0.95) * sensitivity_table_rival_buffer$Cluster.s.e.
sensitivity_table_rival_buffer$Analysis[sensitivity_table_rival_buffer$Analysis=="PA"] <- "No Precip."
sensitivity_table_rival_combined <- read.csv("Outputs/regression_table_total_precip.csv")
sensitivity_table_rival_combined$Lower <- sensitivity_table_rival_combined$Estimate - 1.96 * sensitivity_table_rival_combined$Cluster.s.e.
sensitivity_table_rival_combined$Upper <- sensitivity_table_rival_combined$Estimate + 1.96 * sensitivity_table_rival_combined$Cluster.s.e.
sensitivity_table_rival_combined$Lower90 <- sensitivity_table_rival_combined$Estimate - qnorm(0.95) * sensitivity_table_rival_combined$Cluster.s.e.
sensitivity_table_rival_combined$Upper90 <- sensitivity_table_rival_combined$Estimate + qnorm(0.95) * sensitivity_table_rival_combined$Cluster.s.e.
sensitivity_table_rival_combined$Analysis[sensitivity_table_rival_combined$Analysis=="PA"] <- "No Precip."
sensitivity_table_rival <- sensitivity_table[sensitivity_table$model!="Aggregate",]
sensitivity_table_rival <- sensitivity_table_rival %>% rename(Analysis=model)
sensitivity_table_rival <- rbind(sensitivity_table_rival , sensitivity_table_rival_buffer, sensitivity_table_rival_combined)
sensitivity_table_rival$Location <- factor(sensitivity_table_rival$Location, levels = c("Inside PA", "Buffer 3km", "Inside PA + Buffer 3km"))
sensitivity_table_rival$`Rival Explanation` <- sensitivity_table_rival$Analysis
sensitivity_table_rival2 <- sensitivity_table_rival %>% filter(Analysis !="No Funding, 2004-2011")
rival_plot <- ggplot(data=sensitivity_table_rival2[sensitivity_table_rival2$Variable=="Tourism",], aes(x= Location, y= Estimate))+
geom_hline(aes(yintercept=0), linetype="dashed")+
geom_linerange(aes(ymin = Lower, ymax = Upper,shape=`Rival Explanation` ), position = position_dodge(0.3), color="magenta4")+
geom_linerange(aes(ymin = Lower90, ymax = Upper90,shape=`Rival Explanation`), position = position_dodge(0.3), color="magenta4", linewidth=1.5)+
geom_point(aes(shape=`Rival Explanation`), position = position_dodge(0.3), size=3.5, color="magenta4")+
theme_classic()+
theme(legend.position="top")+
ylab("Estimated Impact of Tourism on Deforestation")
write.csv(sensitivity_table_rival, "Outputs/regression_table_RIVAL.csv")
#90m resolution
sensitivity_resolution <- read.csv("Outputs/points_summary_output.csv")
sensitivity_resolution$Location <- "Inside PA"
sensitivity_resolution_buffer <- read.csv("~/Documents/GitHub/TourismDeforestation/Data/points_summary_output_buffer.csv")
sensitivity_resolution_buffer$Location <- "Buffer 3km"
sensitivity_resolution_combined <- read.csv("~/Documents/GitHub/TourismDeforestation/Data/points_summary_output_combined.csv")
sensitivity_resolution_combined$Location <- "Inside PA + Buffer 3km"
sensitivity_resolution <- rbind(sensitivity_resolution,sensitivity_resolution_buffer,sensitivity_resolution_combined)
sensitivity_resolution$Lower <- sensitivity_resolution$Estimate - 1.96 * sensitivity_resolution$Cluster.s.e.
sensitivity_resolution$Upper <- sensitivity_resolution$Estimate + 1.96 * sensitivity_resolution$Cluster.s.e.
sensitivity_resolution$Lower90 <- sensitivity_resolution$Estimate -qnorm(0.95) * sensitivity_resolution$Cluster.s.e.
sensitivity_resolution$Upper90 <- sensitivity_resolution$Estimate + qnorm(0.95) * sensitivity_resolution$Cluster.s.e.
sensitivity_resolution$Location <- factor(sensitivity_resolution$Location, levels=c("Inside PA", "Buffer 3km", "Inside PA + Buffer 3km"))
resolution_plot <- ggplot(data=sensitivity_resolution[sensitivity_resolution$Variable=="Tourism",], aes(x = Location , y= Estimate))+
geom_hline(aes(yintercept=0), linetype="dashed")+
geom_linerange(aes(ymin = Lower, ymax = Upper), position = position_dodge(0.3),color="magenta4")+
geom_linerange(aes(ymin = Lower90, ymax = Upper90), position = position_dodge(0.3),color="magenta4", linewidth=1.5)+
geom_point(position = position_dodge(0.3), size=3.5, color="magenta4",)+
theme_classic()+
ylab("Estimated Impact of Tourism on Deforestation")
write.csv(sensitivity_resolution, "Outputs/regression_table_RESOLUTION.csv")
#PUT ALL TOGETHER
sensitivity_plots1 <- ggarrange(distance_plot, aggregate_plot, labels = c("a", "b"), nrow=2, ncol=1)
sensitivity_plots2 <- ggarrange(rival_plot, resolution_plot, labels = c("a", "b"), nrow=2, ncol=1)
```