Grauger Causal tests for burned area vs. Water Deficit in 3 Savannas of South Americ"

knitr::opts_chunk$set(fig.width=8, fig.height=4, fig.path = "figure/Granger_")

require(lmtest)
require(ggplot2)

Gran Sabana

load("df_gscru.rda")
#summary(df_gscru)

plot(df_gscru$Date,df_gscru$Deficit,type="b")
plot of chunk unnamed-chunk-2

plot of chunk unnamed-chunk-2

plot(df_gscru$Date,(df_gscru$Surplus + df_gscru$Deficitu)/100,type="b")
lines(df_gscru$Date,df_gscru$BAperc,col="red")
plot of chunk unnamed-chunk-2

plot of chunk unnamed-chunk-2

df_gscru$Deficitperc <- df_gscru$Deficit*100/max(df_gscru$Deficit)
df_gscru$BAperc2 <- df_gscru$BAperc*100/max(df_gscru$BAperc)
plot(df_gscru$Date,df_gscru$Deficitperc,type="b")
lines(df_gscru$Date,df_gscru$BAperc2,col="red")
title("Gran Sabana")
plot of chunk unnamed-chunk-2

plot of chunk unnamed-chunk-2

gtot <- matrix(1:15,ncol=3)
for(ord in 1:5){
  g <- grangertest(x=df_gscru$Deficitperc,
            y=df_gscru$BAperc2, order=ord)
  print(g)
  gtot[ord,1] <- g[2,4]
}
## Granger causality test
## 
## Model 1: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:1) + Lags(df_gscru$Deficitperc, 1:1)
## Model 2: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:1)
##   Res.Df Df      F   Pr(>F)   
## 1    116                      
## 2    117 -1 8.6693 0.003911 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:2) + Lags(df_gscru$Deficitperc, 1:2)
## Model 2: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:2)
##   Res.Df Df      F  Pr(>F)  
## 1    113                    
## 2    115 -2 3.9715 0.02153 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:3) + Lags(df_gscru$Deficitperc, 1:3)
## Model 2: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:3)
##   Res.Df Df      F  Pr(>F)  
## 1    110                    
## 2    113 -3 2.5922 0.05634 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:4) + Lags(df_gscru$Deficitperc, 1:4)
## Model 2: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:4)
##   Res.Df Df      F  Pr(>F)  
## 1    107                    
## 2    111 -4 2.0239 0.09622 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:5) + Lags(df_gscru$Deficitperc, 1:5)
## Model 2: df_gscru$BAperc2 ~ Lags(df_gscru$BAperc2, 1:5)
##   Res.Df Df      F Pr(>F)
## 1    104                 
## 2    109 -5 1.5214 0.1895

Llanos

load("df_llanoscru.rda")
df_llanoscru$Deficitperc <- df_llanoscru$Deficit*100/max(df_llanoscru$Deficit)
df_llanoscru$BAperc2 <- df_llanoscru$BAperc*100/max(df_llanoscru$BAperc)
plot(df_llanoscru$Date,df_llanoscru$Deficitperc,type="b")
lines(df_llanoscru$Date,df_llanoscru$BAperc2,col="red")
title("Llanos")
plot of chunk unnamed-chunk-3

plot of chunk unnamed-chunk-3

for(ord in 1:5){
  g <- grangertest(x=df_llanoscru$Deficitperc,
                   y=df_llanoscru$BAperc2, order=ord)
  print(g)
  gtot[ord,2] <- g[2,4]
}
## Granger causality test
## 
## Model 1: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:1) + Lags(df_llanoscru$Deficitperc, 1:1)
## Model 2: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:1)
##   Res.Df Df      F  Pr(>F)  
## 1    116                    
## 2    117 -1 6.3165 0.01333 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:2) + Lags(df_llanoscru$Deficitperc, 1:2)
## Model 2: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:2)
##   Res.Df Df      F  Pr(>F)  
## 1    113                    
## 2    115 -2 3.2342 0.04306 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:3) + Lags(df_llanoscru$Deficitperc, 1:3)
## Model 2: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:3)
##   Res.Df Df      F  Pr(>F)  
## 1    110                    
## 2    113 -3 3.2976 0.02317 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:4) + Lags(df_llanoscru$Deficitperc, 1:4)
## Model 2: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:4)
##   Res.Df Df     F  Pr(>F)  
## 1    107                   
## 2    111 -4 2.552 0.04321 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:5) + Lags(df_llanoscru$Deficitperc, 1:5)
## Model 2: df_llanoscru$BAperc2 ~ Lags(df_llanoscru$BAperc2, 1:5)
##   Res.Df Df      F  Pr(>F)  
## 1    104                    
## 2    109 -5 1.9962 0.08526 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Beni

load("df_benicru.rda")
df_benicru$Deficitperc <- df_benicru$Deficit*100/max(df_benicru$Deficit)
df_benicru$BAperc2 <- df_benicru$BAperc*100/max(df_benicru$BAperc)
plot(df_benicru$Date,df_benicru$Deficitperc,type="b")
lines(df_benicru$Date,df_benicru$BAperc2,col="red")
title("Beni")
plot of chunk unnamed-chunk-4

plot of chunk unnamed-chunk-4

for(ord in 1:5){
  g <- grangertest(x=df_benicru$Deficitperc,
                   y=df_benicru$BAperc2, order=ord)
  print(g)
  gtot[ord,3] <- g[2,4]
}
## Granger causality test
## 
## Model 1: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:1) + Lags(df_benicru$Deficitperc, 1:1)
## Model 2: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:1)
##   Res.Df Df      F    Pr(>F)    
## 1    116                        
## 2    117 -1 16.563 8.623e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:2) + Lags(df_benicru$Deficitperc, 1:2)
## Model 2: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:2)
##   Res.Df Df      F   Pr(>F)   
## 1    113                      
## 2    115 -2 5.6058 0.004773 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:3) + Lags(df_benicru$Deficitperc, 1:3)
## Model 2: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:3)
##   Res.Df Df      F   Pr(>F)   
## 1    110                      
## 2    113 -3 4.6599 0.004186 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:4) + Lags(df_benicru$Deficitperc, 1:4)
## Model 2: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:4)
##   Res.Df Df      F   Pr(>F)   
## 1    107                      
## 2    111 -4 4.2455 0.003147 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:5) + Lags(df_benicru$Deficitperc, 1:5)
## Model 2: df_benicru$BAperc2 ~ Lags(df_benicru$BAperc2, 1:5)
##   Res.Df Df      F  Pr(>F)   
## 1    104                     
## 2    109 -5 3.4132 0.00677 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Compare among Savannas

matplot(1:5,gtot,xlab="Lag",ylab="p", pch=20,
        main="Significance of Granger Causal Test")
abline(a=0.05, b=0,lty=3)
legend("topleft",legend=c("GS","LL","B"),pch=20,col=1:3)
plot of chunk unnamed-chunk-5

plot of chunk unnamed-chunk-5

#df <- rbind(df_gscru,df_llanoscru,df_benicru)
#Sofia ha puesto Deficitu y Et0 solo en GS !!!
df <- rbind(df_gscru[,-(13:14)],df_llanoscru,df_benicru)
savanna <- c(rep("GS", nrow(df_gscru)), rep("LL", nrow(df_llanoscru)),rep("B", nrow(df_benicru)))
df$savanna <- savanna
df_AQMcru <- df
rm(df)
ggAQM <- ggplot(data=df_AQMcru) +
  geom_line(aes(x=Date,y=Deficitperc)) +
  geom_line(aes(x=Date,y=BAperc2),col="red") +
  ylab("Percent") +
  facet_wrap(~savanna,ncol=1)
print(ggAQM)