Granger Causality Test


##載入套件

library(lmtest)
library(dplyr)
library(ggplot2)
library(openxlsx)
library(vars)
library(tseries)
library(jtools)
library(forecast)
library(urca)


##建立時間序列

season <- read.xlsx("C:/Users/isisg/OneDrive/Desktop/summer program/hdd/new data/byseason1.xlsx", sheet = '內插')

cover.timeseries <- ts(season$cover距平, start = c(1997,1),frequency = 4)
biomass.timeseries <- ts(season$biomass距平, start = c(1997,1),frequency = 4)
oni.timeseries <- ts(season$ONI0, start = c(1997,1),frequency = 4)
pmm.timeseries <- ts(season$PMM0, start = c(1997,1),frequency = 4)
nino4.timeseries <- ts(season$NINO40, start = c(1997,1),frequency = 4)
cp.timeseries <- ts(season$CPENSO0, start = c(1997,1),frequency = 4)
temp.timeseries <- ts(season$TX距平0, start = c(1997,1),frequency = 4)
rain.timeseries <- ts(season$PP距平0, start = c(1997,1),frequency = 4)
salt.timeseries <- ts(season$Salt距平0, start = c(1997,1),frequency = 4)


##檢測平穩性

ADF法

adf.test(na.omit(cover.timeseries))
adf.test(na.omit(biomass.timeseries))

adf.test(na.omit(pmm.timeseries))
adf.test(na.omit(oni.timeseries))
adf.test(na.omit(nino4.timeseries))
adf.test(na.omit(cp.timeseries))

adf.test(na.omit(temp.timeseries))
adf.test(na.omit(rain.timeseries))
adf.test(na.omit(salt.timeseries))
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(cover.timeseries)
## Dickey-Fuller = -1.6921, Lag order = 4, p-value = 0.7022
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(biomass.timeseries)
## Dickey-Fuller = -1.2938, Lag order = 4, p-value = 0.8655
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(pmm.timeseries)
## Dickey-Fuller = -2.5903, Lag order = 4, p-value = 0.3325
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(oni.timeseries)
## Dickey-Fuller = -3.6285, Lag order = 4, p-value = 0.03462
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(nino4.timeseries)
## Dickey-Fuller = -3.4569, Lag order = 4, p-value = 0.04983
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(cp.timeseries)
## Dickey-Fuller = -3.834, Lag order = 4, p-value = 0.02091
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(temp.timeseries)
## Dickey-Fuller = -2.2262, Lag order = 4, p-value = 0.483
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(rain.timeseries)
## Dickey-Fuller = -4.7815, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary
## 
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(salt.timeseries)
## Dickey-Fuller = -2.4348, Lag order = 4, p-value = 0.3975
## alternative hypothesis: stationary


KPSS法

kpss.test(na.omit(cover.timeseries))
kpss.test(na.omit(biomass.timeseries))

kpss.test(na.omit(pmm.timeseries))
kpss.test(na.omit(oni.timeseries))
kpss.test(na.omit(nino4.timeseries))
kpss.test(na.omit(cp.timeseries))

kpss.test(na.omit(temp.timeseries))
kpss.test(na.omit(rain.timeseries))
kpss.test(na.omit(salt.timeseries))
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(cover.timeseries)
## KPSS Level = 0.8314, Truncation lag parameter = 3, p-value = 0.01
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(biomass.timeseries)
## KPSS Level = 0.55658, Truncation lag parameter = 3, p-value = 0.02892
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(pmm.timeseries)
## KPSS Level = 0.8593, Truncation lag parameter = 3, p-value = 0.01
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(oni.timeseries)
## KPSS Level = 0.11392, Truncation lag parameter = 3, p-value = 0.1
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(nino4.timeseries)
## KPSS Level = 0.19939, Truncation lag parameter = 3, p-value = 0.1
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(cp.timeseries)
## KPSS Level = 0.22521, Truncation lag parameter = 3, p-value = 0.1
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(temp.timeseries)
## KPSS Level = 0.78001, Truncation lag parameter = 3, p-value = 0.01
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(rain.timeseries)
## KPSS Level = 0.12981, Truncation lag parameter = 3, p-value = 0.1
## 
## 
##  KPSS Test for Level Stationarity
## 
## data:  na.omit(salt.timeseries)
## KPSS Level = 0.17933, Truncation lag parameter = 3, p-value = 0.1


##檢測差分

ADF法

ndiffs(cover.timeseries, test="adf")
ndiffs(biomass.timeseries, test="adf")

ndiffs(pmm.timeseries, test="adf")
ndiffs(oni.timeseries, test="adf")
ndiffs(nino4.timeseries, test="adf")
ndiffs(cp.timeseries, test="adf")

ndiffs(temp.timeseries, test="adf")
ndiffs(rain.timeseries, test="adf")
ndiffs(salt.timeseries, test="adf")
## [1] 1
## [1] 1
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0


KPSS法

ndiffs(cover.timeseries, test="kpss")
ndiffs(biomass.timeseries, test="kpss")

ndiffs(pmm.timeseries, test="kpss")
ndiffs(oni.timeseries, test="kpss")
ndiffs(nino4.timeseries, test="kpss")
ndiffs(cp.timeseries, test="kpss")

ndiffs(temp.timeseries, test="kpss")
ndiffs(rain.timeseries, test="kpss")
ndiffs(salt.timeseries, test="kpss")
## [1] 1
## [1] 1
## [1] 1
## [1] 0
## [1] 0
## [1] 0
## [1] 1
## [1] 0
## [1] 0

cover, biomass必做差分;PMM, temp可做可不做

##檢測季節性差分

nsdiffs(cover.timeseries)
nsdiffs(biomass.timeseries)

nsdiffs(pmm.timeseries)
nsdiffs(oni.timeseries)
nsdiffs(nino4.timeseries)
nsdiffs(cp.timeseries)

nsdiffs(temp.timeseries)
nsdiffs(rain.timeseries)
nsdiffs(salt.timeseries)
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0
## [1] 0



##差分比較

coverdiff <- diff(cover.timeseries)

par(mfrow=c(2,2))
plot.ts(cover.timeseries,main="原始cover")
plot.ts(coverdiff,main="一階差分cover")
acf(cover.timeseries, na = na.omit, ,main="原始cover"); acf(coverdiff, na = na.omit, main="一階差分cover")

biomassdiff <- diff(biomass.timeseries)
par(mfrow=c(2,2))
plot.ts(biomass.timeseries,main="原始biomass")
plot.ts(coverdiff,main="一階差分biomass")
acf(biomass.timeseries, na = na.omit, ,main="原始biomass");acf(coverdiff, na = na.omit, main="一階差分biomass")

pmmdiff <- diff(pmm.timeseries)
#par(mfrow=c(2,2))
#plot.ts(pmm.timeseries,main="原始pmm")
#plot.ts(pmmdiff,main="一階差分pmm")
#acf(pmm.timeseries, na = na.omit, ,main="原始pmm");acf(pmmdiff, na = na.omit, main="一階差分pmm")
tempdiff <- diff(temp.timeseries)

#par(mfrow=c(2,2))
#plot.ts(temp.timeseries,main="原始temp")
#plot.ts(tempdiff,main="一階差分temp")
#(temp.timeseries, na = na.omit, ,main="原始temp");acf(tempdiff, na = na.omit, main="一階差分temp")



##建立校準時間序列

因為差分的原理為第一項與第二項之差為新序列之第一項,因此計算a與b之格蘭傑因果時,其結果會與現實差了一個時間單位。然而當此差距會包含0時,將導致一盲區無法被計算(沒有自回歸),因此需要建立一個向後推移一單位之時間序列。

cover1.timeseries <- ts(season$cover距平1, start = c(1997,1),frequency = 4)
biomass1.timeseries <- ts(season$biomass距平1, start = c(1997,1),frequency = 4)
pmm1.timeseries <- ts(season$PMM1, start = c(1997,1),frequency = 4)
temp1.timeseries <- ts(season$TX距平1, start = c(1997,1),frequency = 4)



##建立校準差分序列

cover1diff <- diff(cover1.timeseries)
biomass1diff <- diff(biomass1.timeseries)
pmm1diff <- diff(pmm1.timeseries)
temp1diff <- diff(temp1.timeseries)



##Biomass與PMM的滯後選擇 & Granger Causality Test

BIOMASS差分與PMM

bio_pmm <- cbind(pmm.timeseries, biomassdiff)
VARselect(y=na.omit(bio_pmm), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      5      1      1      5 
## 
## $criteria
##                 1          2          3          4          5          6
## AIC(n) -1.0142979 -1.0510566 -1.0587417 -1.0573471 -1.0670205 -0.9913257
## HQ(n)  -0.9402701 -0.9276770 -0.8860102 -0.8352637 -0.7955852 -0.6705385
## SC(n)  -0.8288989 -0.7420582 -0.6261439 -0.5011500 -0.3872240 -0.1879298
## FPE(n)  0.3626879  0.3497067  0.3472703  0.3481844  0.3455014  0.3737179
##                  7          8
## AIC(n) -0.90715411 -0.8607084
## HQ(n)  -0.53701509 -0.4412175
## SC(n)   0.01984113  0.1898862
## FPE(n)  0.40810781  0.4296932

AIC best at 5, BIC(SC) best at 1.

grangertest(biomassdiff ~ pmm.timeseries, order = 1, data = bio_pmm)
grangertest(pmm.timeseries ~ biomassdiff, order = 1, data = bio_pmm)
grangertest(biomassdiff ~ pmm.timeseries, order = 2, data = bio_pmm)
grangertest(pmm.timeseries ~ biomassdiff, order = 2, data = bio_pmm)
grangertest(biomassdiff ~ pmm.timeseries, order = 3, data = bio_pmm)
grangertest(pmm.timeseries ~ biomassdiff, order = 3, data = bio_pmm)
grangertest(biomassdiff ~ pmm.timeseries, order = 4, data = bio_pmm)
grangertest(pmm.timeseries ~ biomassdiff, order = 4, data = bio_pmm)
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:1) + Lags(pmm.timeseries, 1:1)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:1)
##   Res.Df Df     F Pr(>F)
## 1     79                
## 2     80 -1 0.164 0.6866
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:1) + Lags(biomassdiff, 1:1)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)  
## 1     79                   
## 2     80 -1 3.8065 0.0546 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:2) + Lags(pmm.timeseries, 1:2)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.0236 0.9767
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:2) + Lags(biomassdiff, 1:2)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 1.7446 0.1816
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:3) + Lags(pmm.timeseries, 1:3)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.1011 0.9592
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:3) + Lags(biomassdiff, 1:3)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.8647 0.4634
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:4) + Lags(pmm.timeseries, 1:4)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.6346 0.6395
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:4) + Lags(biomassdiff, 1:4)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.8408 0.1307

Biomass lead PMM by 1. (實質上Biomass與領先1.5季)

BIOMASS1差分與PMM

bio1_pmm <- cbind(pmm.timeseries, biomass1diff)
VARselect(y=na.omit(bio1_pmm), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      3      2      1      3 
## 
## $criteria
##                 1          2          3          4          5          6
## AIC(n) -0.9235428 -1.0423719 -1.0881076 -1.0168086 -0.9951740 -0.9836241
## HQ(n)  -0.8495150 -0.9189922 -0.9153761 -0.7947252 -0.7237388 -0.6628370
## SC(n)  -0.7381438 -0.7333735 -0.6555098 -0.4606114 -0.3153775 -0.1802282
## FPE(n)  0.3971436  0.3527570  0.3372206  0.3625893  0.3712379  0.3766072
##                  7          8
## AIC(n) -0.90363397 -0.8821845
## HQ(n)  -0.53349495 -0.4626936
## SC(n)   0.02336128  0.1684101
## FPE(n)  0.40954694  0.4205634

AIC best at 3, BIC(SC) best at 1.

grangertest(biomass1diff ~ pmm.timeseries, order = 1, data = bio1_pmm)
grangertest(pmm.timeseries ~ biomass1diff, order = 1, data = bio1_pmm)
grangertest(biomass1diff ~ pmm.timeseries, order = 2, data = bio1_pmm)
grangertest(pmm.timeseries ~ biomass1diff, order = 2, data = bio1_pmm)
grangertest(biomass1diff ~ pmm.timeseries, order = 3, data = bio1_pmm)
grangertest(pmm.timeseries ~ biomass1diff, order = 3, data = bio1_pmm)
grangertest(biomass1diff ~ pmm.timeseries, order = 4, data = bio1_pmm)
grangertest(pmm.timeseries ~ biomass1diff, order = 4, data = bio1_pmm)
grangertest(biomass1diff ~ pmm.timeseries, order = 5, data = bio1_pmm)
grangertest(pmm.timeseries ~ biomass1diff, order = 5, data = bio1_pmm)
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:1) + Lags(pmm.timeseries, 1:1)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:1)
##   Res.Df Df  F Pr(>F)
## 1     79             
## 2     80 -1  0 0.9959
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:1) + Lags(biomass1diff, 1:1)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 2.3891 0.1262
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:2) + Lags(pmm.timeseries, 1:2)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:2)
##   Res.Df Df     F Pr(>F)
## 1     76                
## 2     78 -2 0.055 0.9465
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:2) + Lags(biomass1diff, 1:2)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:2)
##   Res.Df Df      F  Pr(>F)  
## 1     76                    
## 2     78 -2 3.9938 0.02243 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:3) + Lags(pmm.timeseries, 1:3)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.8312  0.481
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:3) + Lags(biomass1diff, 1:3)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:3)
##   Res.Df Df      F  Pr(>F)  
## 1     73                    
## 2     76 -3 3.2092 0.02793 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:4) + Lags(pmm.timeseries, 1:4)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.6521 0.6273
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:4) + Lags(biomass1diff, 1:4)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:4)
##   Res.Df Df      F  Pr(>F)  
## 1     70                    
## 2     74 -4 2.3539 0.06216 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:5) + Lags(pmm.timeseries, 1:5)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 0.6801 0.6401
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:5) + Lags(biomass1diff, 1:5)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:5)
##   Res.Df Df      F  Pr(>F)  
## 1     67                    
## 2     72 -5 2.3023 0.05432 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Biomass lead PMM by 3. (實質上Biomass領先PMM2季)

##Biomass與Temp的滯後選擇 & Granger

Biomass差分與Temp

temp_biomass <- cbind(temp.timeseries, biomassdiff)
VARselect(y=na.omit(temp_biomass), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      5      1      1      5 
## 
## $criteria
##                  1           2           3           4           5           6
## AIC(n) -3.64030637 -3.66062270 -3.66806190 -3.67286969 -3.69050411 -3.65763917
## HQ(n)  -3.56627857 -3.53724303 -3.49533036 -3.45078627 -3.41906882 -3.33685202
## SC(n)  -3.45490732 -3.35162428 -3.23546412 -3.11667254 -3.01070759 -2.85424329
## FPE(n)  0.02624654  0.02572668  0.02555373  0.02546257  0.02506601  0.02597638
##                  7           8
## AIC(n) -3.56761590 -3.49338562
## HQ(n)  -3.19747688 -3.07389473
## SC(n)  -2.64062066 -2.44279101
## FPE(n)  0.02853324  0.03088881

AIC best at 5, BIC best at 1.

grangertest(biomassdiff ~ temp.timeseries, order = 1, data = temp_biomass)
grangertest(temp.timeseries ~ biomassdiff, order = 1, data = temp_biomass)
grangertest(biomassdiff ~ temp.timeseries, order = 2, data = temp_biomass)
grangertest(temp.timeseries ~ biomassdiff, order = 2, data = temp_biomass)
grangertest(biomassdiff ~ temp.timeseries, order = 3, data = temp_biomass)
grangertest(temp.timeseries ~ biomassdiff, order = 3, data = temp_biomass)
grangertest(biomassdiff ~ temp.timeseries, order = 4, data = temp_biomass)
grangertest(temp.timeseries ~ biomassdiff, order = 4, data = temp_biomass)
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:1) + Lags(temp.timeseries, 1:1)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.0149 0.9031
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:1) + Lags(biomassdiff, 1:1)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:1)
##   Res.Df Df     F Pr(>F)
## 1     79                
## 2     80 -1 0.112 0.7388
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:2) + Lags(temp.timeseries, 1:2)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.4357 0.6484
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:2) + Lags(biomassdiff, 1:2)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.6661 0.5167
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:3) + Lags(temp.timeseries, 1:3)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.3548 0.7858
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:3) + Lags(biomassdiff, 1:3)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.7378 0.5329
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:4) + Lags(temp.timeseries, 1:4)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:4)
##   Res.Df Df     F Pr(>F)
## 1     70                
## 2     74 -4 0.307 0.8724
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:4) + Lags(biomassdiff, 1:4)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.6976 0.1603



Biomass1差分與Temp

temp_biomass1 <- cbind(temp.timeseries, biomass1diff)
VARselect(y=na.omit(temp_biomass1), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      3      1      1      3 
## 
## $criteria
##                  1           2           3           4           5           6
## AIC(n) -3.64759207 -3.65713865 -3.67318987 -3.60843188 -3.66587822 -3.67300128
## HQ(n)  -3.57356426 -3.53375897 -3.50045833 -3.38634847 -3.39444293 -3.35221412
## SC(n)  -3.46219302 -3.34814023 -3.24059209 -3.05223473 -2.98608170 -2.86960540
## FPE(n)  0.02605601  0.02581647  0.02542302  0.02715734  0.02569095  0.02558038
##                  7           8
## AIC(n) -3.59927194 -3.50179199
## HQ(n)  -3.22913292 -3.08230109
## SC(n)  -2.67227670 -2.45119737
## FPE(n)  0.02764413  0.03063024

AIC best at 5, BIC best at 1.

grangertest(biomass1diff ~ temp.timeseries, order = 1, data = temp_biomass1)
grangertest(temp.timeseries ~ biomass1diff, order = 1, data = temp_biomass1)
grangertest(biomass1diff ~ temp.timeseries, order = 2, data = temp_biomass1)
grangertest(temp.timeseries ~ biomass1diff, order = 2, data = temp_biomass1)
grangertest(biomass1diff ~ temp.timeseries, order = 3, data = temp_biomass1)
grangertest(temp.timeseries ~ biomass1diff, order = 3, data = temp_biomass1)
grangertest(biomass1diff ~ temp.timeseries, order = 4, data = temp_biomass1)
grangertest(temp.timeseries ~ biomass1diff, order = 4, data = temp_biomass1)
grangertest(biomass1diff ~ temp.timeseries, order = 5, data = temp_biomass1)
grangertest(temp.timeseries ~ biomass1diff, order = 5, data = temp_biomass1)
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:1) + Lags(temp.timeseries, 1:1)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.4509 0.5039
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:1) + Lags(biomass1diff, 1:1)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.3418 0.5604
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:2) + Lags(temp.timeseries, 1:2)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.2952 0.7452
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:2) + Lags(biomass1diff, 1:2)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.4851 0.6175
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:3) + Lags(temp.timeseries, 1:3)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.2271 0.8773
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:3) + Lags(biomass1diff, 1:3)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.6559 0.5818
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:4) + Lags(temp.timeseries, 1:4)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.1745 0.9508
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:4) + Lags(biomass1diff, 1:4)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.6836 0.6057
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:5) + Lags(temp.timeseries, 1:5)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 0.2487 0.9391
## Granger causality test
## 
## Model 1: temp.timeseries ~ Lags(temp.timeseries, 1:5) + Lags(biomass1diff, 1:5)
## Model 2: temp.timeseries ~ Lags(temp.timeseries, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 1.1218 0.3573

No significant result.

##Biomass與Salt的滯後選擇 & Granger Causality Test

BIOMASS差分與Salt

bio_salt <- cbind(salt.timeseries, biomassdiff)
VARselect(y=na.omit(bio_salt), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      2      2      1      2 
## 
## $criteria
##                 1          2          3          4          5          6
## AIC(n) -2.0013017 -2.1037585 -2.0286503 -2.0782186 -2.0553771 -2.0522108
## HQ(n)  -1.9272739 -1.9803788 -1.8559188 -1.8561352 -1.7839418 -1.7314236
## SC(n)  -1.8159026 -1.7947601 -1.5960526 -1.5220214 -1.3755806 -1.2488149
## FPE(n)  0.1351708  0.1220454  0.1316563  0.1254442  0.1285914  0.1293622
##                 7          8
## AIC(n) -2.0603694 -2.0663946
## HQ(n)  -1.6902304 -1.6469037
## SC(n)  -1.1333742 -1.0158000
## FPE(n)  0.1288071  0.1286873

AIC, BIC(SC) best at 2.

grangertest(biomassdiff ~ salt.timeseries, order = 1, data = bio_salt)
grangertest(salt.timeseries ~ biomassdiff, order = 1, data = bio_salt)
grangertest(biomassdiff ~ salt.timeseries, order = 2, data = bio_salt)
grangertest(salt.timeseries ~ biomassdiff, order = 2, data = bio_salt)
grangertest(biomassdiff ~ salt.timeseries, order = 3, data = bio_salt)
grangertest(salt.timeseries ~ biomassdiff, order = 3, data = bio_salt)
grangertest(biomassdiff ~ salt.timeseries, order = 4, data = bio_salt)
grangertest(salt.timeseries ~ biomassdiff, order = 4, data = bio_salt)
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:1) + Lags(salt.timeseries, 1:1)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.9089 0.3433
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:1) + Lags(biomassdiff, 1:1)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.0019 0.9653
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:2) + Lags(salt.timeseries, 1:2)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.8818 0.4182
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:2) + Lags(biomassdiff, 1:2)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:2)
##   Res.Df Df     F Pr(>F)
## 1     76                
## 2     78 -2 0.931 0.3986
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:3) + Lags(salt.timeseries, 1:3)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.7218 0.5422
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:3) + Lags(biomassdiff, 1:3)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.8013 0.4972
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:4) + Lags(salt.timeseries, 1:4)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.8516 0.4974
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:4) + Lags(biomassdiff, 1:4)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.5973 0.6658

No significant result.

BIOMASS1差分與Salt

bio1_salt <- cbind(salt.timeseries, biomass1diff)
VARselect(y=na.omit(bio1_salt), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      2      2      1      2 
## 
## $criteria
##                 1          2          3          4          5          6
## AIC(n) -1.9832825 -2.0466255 -1.9994639 -1.9809571 -1.9258896 -1.8492010
## HQ(n)  -1.9092547 -1.9232458 -1.8267324 -1.7588737 -1.6544543 -1.5284138
## SC(n)  -1.7978835 -1.7376271 -1.5668661 -1.4247599 -1.2460931 -1.0458051
## FPE(n)  0.1376285  0.1292212  0.1355555  0.1382581  0.1463685  0.1584796
##                 7          8
## AIC(n) -1.8517057 -1.8262288
## HQ(n)  -1.4815667 -1.4067379
## SC(n)  -0.9247105 -0.7756342
## FPE(n)  0.1586943  0.1636207

AIC best at 1, BIC(SC) best at 2.

grangertest(biomass1diff ~ salt.timeseries, order = 1, data = bio1_salt)
grangertest(salt.timeseries ~ biomass1diff, order = 1, data = bio1_salt)
grangertest(biomass1diff ~ salt.timeseries, order = 2, data = bio1_salt)
grangertest(salt.timeseries ~ biomass1diff, order = 2, data = bio1_salt)
grangertest(biomass1diff ~ salt.timeseries, order = 3, data = bio1_salt)
grangertest(salt.timeseries ~ biomass1diff, order = 3, data = bio1_salt)
grangertest(biomass1diff ~ salt.timeseries, order = 4, data = bio1_salt)
grangertest(salt.timeseries ~ biomass1diff, order = 4, data = bio1_salt)
grangertest(biomass1diff ~ salt.timeseries, order = 5, data = bio1_salt)
grangertest(salt.timeseries ~ biomass1diff, order = 5, data = bio1_salt)
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:1) + Lags(salt.timeseries, 1:1)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.8825 0.3504
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:1) + Lags(biomass1diff, 1:1)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.0548 0.8155
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:2) + Lags(salt.timeseries, 1:2)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.8387 0.4363
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:2) + Lags(biomass1diff, 1:2)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.3051 0.7379
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:3) + Lags(salt.timeseries, 1:3)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.6707 0.5727
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:3) + Lags(biomass1diff, 1:3)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:3)
##   Res.Df Df     F Pr(>F)
## 1     73                
## 2     76 -3 0.559 0.6437
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:4) + Lags(salt.timeseries, 1:4)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.9022 0.4675
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:4) + Lags(biomass1diff, 1:4)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:4)
##   Res.Df Df     F Pr(>F)
## 1     70                
## 2     74 -4 0.409 0.8015
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:5) + Lags(salt.timeseries, 1:5)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 0.9464  0.457
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:5) + Lags(biomass1diff, 1:5)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 0.4916 0.7814

No significant result.

BIOMASS差分與Rain

bio_rain <- cbind(rain.timeseries, biomassdiff)
VARselect(y=na.omit(bio_rain), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      2      2      1      2 
## 
## $criteria
##                 1          2          3          4          5          6
## AIC(n)   6.701146   6.628116   6.694815   6.738035   6.725236   6.823009
## HQ(n)    6.775174   6.751495   6.867546   6.960118   6.996672   7.143796
## SC(n)    6.886545   6.937114   7.127412   7.294232   7.405033   7.626405
## FPE(n) 813.406914 756.355750 809.085320 845.864264 836.727150 925.265050
##                 7           8
## AIC(n)   6.873256    6.909575
## HQ(n)    7.243395    7.329066
## SC(n)    7.800251    7.960170
## FPE(n) 976.706544 1018.004288

AIC, BIC(SC) best at 2.

grangertest(biomassdiff ~ rain.timeseries, order = 1, data = bio_rain)
grangertest(rain.timeseries ~ biomassdiff, order = 1, data = bio_rain)
grangertest(biomassdiff ~ rain.timeseries, order = 2, data = bio_rain)
grangertest(rain.timeseries ~ biomassdiff, order = 2, data = bio_rain)
grangertest(biomassdiff ~ rain.timeseries, order = 3, data = bio_rain)
grangertest(rain.timeseries ~ biomassdiff, order = 3, data = bio_rain)
grangertest(biomassdiff ~ rain.timeseries, order = 4, data = bio_rain)
grangertest(rain.timeseries ~ biomassdiff, order = 4, data = bio_rain)
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:1) + Lags(rain.timeseries, 1:1)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.4933 0.4845
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:1) + Lags(biomassdiff, 1:1)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.7246 0.3972
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:2) + Lags(rain.timeseries, 1:2)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:2)
##   Res.Df Df      F  Pr(>F)  
## 1     76                    
## 2     78 -2 2.5968 0.08112 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:2) + Lags(biomassdiff, 1:2)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 1.9967 0.1428
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:3) + Lags(rain.timeseries, 1:3)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.7129 0.1718
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:3) + Lags(biomassdiff, 1:3)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.9364 0.1312
## Granger causality test
## 
## Model 1: biomassdiff ~ Lags(biomassdiff, 1:4) + Lags(rain.timeseries, 1:4)
## Model 2: biomassdiff ~ Lags(biomassdiff, 1:4)
##   Res.Df Df     F Pr(>F)
## 1     70                
## 2     74 -4 1.946 0.1124
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:4) + Lags(biomassdiff, 1:4)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.4631 0.2227



BIOMASS1差分與rain

bio1_rain <- cbind(rain.timeseries, biomass1diff)
VARselect(y=na.omit(bio1_rain), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      1      1      1      1 
## 
## $criteria
##                 1          2          3          4          5          6
## AIC(n)   6.654225   6.656723   6.722605   6.736654   6.783025   6.829551
## HQ(n)    6.728253   6.780103   6.895337   6.958737   7.054460   7.150338
## SC(n)    6.839624   6.965721   7.155203   7.292851   7.462822   7.632947
## FPE(n) 776.122619 778.305553 831.885570 844.696799 886.505144 931.338478
##                 7           8
## AIC(n)   6.831679    6.925184
## HQ(n)    7.201818    7.344674
## SC(n)    7.758674    7.975778
## FPE(n) 936.930541 1034.018560

AIC best at 1, BIC(SC) best at 2.

grangertest(biomass1diff ~ rain.timeseries, order = 1, data = bio1_rain)
grangertest(rain.timeseries ~ biomass1diff, order = 1, data = bio1_rain)
grangertest(biomass1diff ~ rain.timeseries, order = 2, data = bio1_rain)
grangertest(rain.timeseries ~ biomass1diff, order = 2, data = bio1_rain)
grangertest(biomass1diff ~ rain.timeseries, order = 3, data = bio1_rain)
grangertest(rain.timeseries ~ biomass1diff, order = 3, data = bio1_rain)
grangertest(biomass1diff ~ rain.timeseries, order = 4, data = bio1_rain)
grangertest(rain.timeseries ~ biomass1diff, order = 4, data = bio1_rain)
grangertest(biomass1diff ~ rain.timeseries, order = 5, data = bio1_rain)
grangertest(rain.timeseries ~ biomass1diff, order = 5, data = bio1_rain)
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:1) + Lags(rain.timeseries, 1:1)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:1)
##   Res.Df Df      F  Pr(>F)  
## 1     79                    
## 2     80 -1 5.7596 0.01875 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:1) + Lags(biomass1diff, 1:1)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 1.2499  0.267
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:2) + Lags(rain.timeseries, 1:2)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 2.3424  0.103
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:2) + Lags(biomass1diff, 1:2)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 1.3759 0.2588
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:3) + Lags(rain.timeseries, 1:3)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.5298  0.214
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:3) + Lags(biomass1diff, 1:3)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.5804 0.2014
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:4) + Lags(rain.timeseries, 1:4)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.6006 0.1838
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:4) + Lags(biomass1diff, 1:4)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.6041 0.1829
## Granger causality test
## 
## Model 1: biomass1diff ~ Lags(biomass1diff, 1:5) + Lags(rain.timeseries, 1:5)
## Model 2: biomass1diff ~ Lags(biomass1diff, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 1.5772 0.1785
## Granger causality test
## 
## Model 1: rain.timeseries ~ Lags(rain.timeseries, 1:5) + Lags(biomass1diff, 1:5)
## Model 2: rain.timeseries ~ Lags(rain.timeseries, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 1.3723 0.2459

##cover與PMM的滯後選擇 & Granger Causality Test

cover差分與PMM

bio_pmm <- cbind(pmm.timeseries, coverdiff)
VARselect(y=na.omit(bio_pmm), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      5      3      1      5 
## 
## $criteria
##                1         2         3         4         5         6         7
## AIC(n)  4.337327  4.346972  4.226343  4.242512  4.209363  4.292804  4.340562
## HQ(n)   4.411355  4.470352  4.399074  4.464595  4.480798  4.613591  4.710701
## SC(n)   4.522726  4.655971  4.658940  4.798709  4.889160  5.096199  5.267558
## FPE(n) 76.509288 77.274835 68.541017 69.744300 67.601135 73.690592 77.594191
##                8
## AIC(n)  4.391809
## HQ(n)   4.811300
## SC(n)   5.442404
## FPE(n) 82.091429

AIC best at 5, BIC(SC) best at 1.

grangertest(coverdiff ~ pmm.timeseries, order = 1, data = bio_pmm)
grangertest(pmm.timeseries ~ coverdiff, order = 1, data = bio_pmm)
grangertest(coverdiff ~ pmm.timeseries, order = 2, data = bio_pmm)
grangertest(pmm.timeseries ~ coverdiff, order = 2, data = bio_pmm)
grangertest(coverdiff ~ pmm.timeseries, order = 3, data = bio_pmm)
grangertest(pmm.timeseries ~ coverdiff, order = 3, data = bio_pmm)
grangertest(coverdiff ~ pmm.timeseries, order = 4, data = bio_pmm)
grangertest(pmm.timeseries ~ coverdiff, order = 4, data = bio_pmm)
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:1) + Lags(pmm.timeseries, 1:1)
## Model 2: coverdiff ~ Lags(coverdiff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.0188 0.8913
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:1) + Lags(coverdiff, 1:1)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:1)
##   Res.Df Df  F Pr(>F)
## 1     79             
## 2     80 -1  0 0.9998
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:2) + Lags(pmm.timeseries, 1:2)
## Model 2: coverdiff ~ Lags(coverdiff, 1:2)
##   Res.Df Df     F Pr(>F)
## 1     76                
## 2     78 -2 0.354  0.703
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:2) + Lags(coverdiff, 1:2)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.2292 0.7957
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:3) + Lags(pmm.timeseries, 1:3)
## Model 2: coverdiff ~ Lags(coverdiff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.9413 0.4252
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:3) + Lags(coverdiff, 1:3)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.5152 0.2177
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:4) + Lags(pmm.timeseries, 1:4)
## Model 2: coverdiff ~ Lags(coverdiff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.8514 0.4975
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:4) + Lags(coverdiff, 1:4)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.0845 0.3709

cover lead PMM by 1 and 5. (實質上cover與PMM同時/領先4季)

cover1差分與PMM

bio1_pmm <- cbind(pmm.timeseries, cover1diff)
VARselect(y=na.omit(bio1_pmm), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      4      1      1      4 
## 
## $criteria
##                1         2         3         4         5         6         7
## AIC(n)  4.329266  4.342998  4.250210  4.248838  4.311245  4.304390  4.375736
## HQ(n)   4.403293  4.466378  4.422941  4.470922  4.582680  4.625177  4.745875
## SC(n)   4.514665  4.651996  4.682807  4.805036  4.991042  5.107786  5.302731
## FPE(n) 75.895018 76.968340 70.196567 70.186961 74.851557 74.549391 80.371998
##                8
## AIC(n)  4.405685
## HQ(n)   4.825176
## SC(n)   5.456280
## FPE(n) 83.238458

AIC best at 4, BIC(SC) best at 1.

grangertest(cover1diff ~ pmm.timeseries, order = 1, data = bio1_pmm)
grangertest(pmm.timeseries ~ cover1diff, order = 1, data = bio1_pmm)
grangertest(cover1diff ~ pmm.timeseries, order = 2, data = bio1_pmm)
grangertest(pmm.timeseries ~ cover1diff, order = 2, data = bio1_pmm)
grangertest(cover1diff ~ pmm.timeseries, order = 3, data = bio1_pmm)
grangertest(pmm.timeseries ~ cover1diff, order = 3, data = bio1_pmm)
grangertest(cover1diff ~ pmm.timeseries, order = 4, data = bio1_pmm)
grangertest(pmm.timeseries ~ cover1diff, order = 4, data = bio1_pmm)
grangertest(cover1diff ~ pmm.timeseries, order = 5, data = bio1_pmm)
grangertest(pmm.timeseries ~ cover1diff, order = 5, data = bio1_pmm)
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:1) + Lags(pmm.timeseries, 1:1)
## Model 2: cover1diff ~ Lags(cover1diff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.4728 0.4937
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:1) + Lags(cover1diff, 1:1)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:1)
##   Res.Df Df      F  Pr(>F)  
## 1     79                    
## 2     80 -1 2.8061 0.09786 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:2) + Lags(pmm.timeseries, 1:2)
## Model 2: cover1diff ~ Lags(cover1diff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.6184 0.5415
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:2) + Lags(cover1diff, 1:2)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 1.9906 0.1437
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:3) + Lags(pmm.timeseries, 1:3)
## Model 2: cover1diff ~ Lags(cover1diff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 0.5372 0.6582
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:3) + Lags(cover1diff, 1:3)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.9201 0.1338
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:4) + Lags(pmm.timeseries, 1:4)
## Model 2: cover1diff ~ Lags(cover1diff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.5734 0.6828
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:4) + Lags(cover1diff, 1:4)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.9923 0.1051
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:5) + Lags(pmm.timeseries, 1:5)
## Model 2: cover1diff ~ Lags(cover1diff, 1:5)
##   Res.Df Df     F Pr(>F)
## 1     67                
## 2     72 -5 0.552 0.7362
## Granger causality test
## 
## Model 1: pmm.timeseries ~ Lags(pmm.timeseries, 1:5) + Lags(cover1diff, 1:5)
## Model 2: pmm.timeseries ~ Lags(pmm.timeseries, 1:5)
##   Res.Df Df     F Pr(>F)
## 1     67                
## 2     72 -5 1.622 0.1662

cover lead PMM by 1. (實質上cover與PMM同時)

##cover與Salt的滯後選擇 & Granger Causality Test

cover差分與Salt

bio_salt <- cbind(salt.timeseries, coverdiff)
VARselect(y=na.omit(bio_salt), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      4      1      1      4 
## 
## $criteria
##                1         2         3         4         5         6         7
## AIC(n)  3.337869  3.299136  3.263972  3.239438  3.296094  3.358735  3.386687
## HQ(n)   3.411896  3.422515  3.436704  3.461522  3.567529  3.679522  3.756826
## SC(n)   3.523268  3.608134  3.696570  3.795635  3.975890  4.162130  4.313682
## FPE(n) 28.161450 27.099951 26.181732 25.578764 27.122277 28.956829 29.892772
##                8
## AIC(n)  3.429431
## HQ(n)   3.848922
## SC(n)   4.480025
## FPE(n) 31.357551

AIC, BIC(SC) best at 2.

grangertest(coverdiff ~ salt.timeseries, order = 1, data = bio_salt)
grangertest(salt.timeseries ~ coverdiff, order = 1, data = bio_salt)
grangertest(coverdiff ~ salt.timeseries, order = 2, data = bio_salt)
grangertest(salt.timeseries ~ coverdiff, order = 2, data = bio_salt)
grangertest(coverdiff ~ salt.timeseries, order = 3, data = bio_salt)
grangertest(salt.timeseries ~ coverdiff, order = 3, data = bio_salt)
grangertest(coverdiff ~ salt.timeseries, order = 4, data = bio_salt)
grangertest(salt.timeseries ~ coverdiff, order = 4, data = bio_salt)
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:1) + Lags(salt.timeseries, 1:1)
## Model 2: coverdiff ~ Lags(coverdiff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.3555 0.5527
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:1) + Lags(coverdiff, 1:1)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:1)
##   Res.Df Df     F Pr(>F)
## 1     79                
## 2     80 -1 0.034 0.8542
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:2) + Lags(salt.timeseries, 1:2)
## Model 2: coverdiff ~ Lags(coverdiff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.0579 0.9438
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:2) + Lags(coverdiff, 1:2)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 1.4071 0.2512
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:3) + Lags(salt.timeseries, 1:3)
## Model 2: coverdiff ~ Lags(coverdiff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.0273 0.3856
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:3) + Lags(coverdiff, 1:3)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.3846 0.2543
## Granger causality test
## 
## Model 1: coverdiff ~ Lags(coverdiff, 1:4) + Lags(salt.timeseries, 1:4)
## Model 2: coverdiff ~ Lags(coverdiff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.1435 0.3434
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:4) + Lags(coverdiff, 1:4)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.5167 0.7236

No significant result.

cover1差分與Salt

bio1_salt <- cbind(salt.timeseries, cover1diff)
VARselect(y=na.omit(bio1_salt), lag.max = 8, type = c("const"))
## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      3      1      1      3 
## 
## $criteria
##                1         2         3         4         5         6         7
## AIC(n)  3.379094  3.332888  3.290836  3.342437  3.411806  3.446901  3.485658
## HQ(n)   3.453122  3.456268  3.463568  3.564521  3.683241  3.767688  3.855797
## SC(n)   3.564493  3.641887  3.723434  3.898635  4.091602  4.250297  4.412653
## FPE(n) 29.346687 28.030253 26.894611 28.353812 30.449437 31.625786 33.002668
##                8
## AIC(n)  3.495164
## HQ(n)   3.914655
## SC(n)   4.545759
## FPE(n) 33.488041

AIC best at 1, BIC(SC) best at 3.

grangertest(cover1diff ~ salt.timeseries, order = 1, data = bio1_salt)
grangertest(salt.timeseries ~ cover1diff, order = 1, data = bio1_salt)
grangertest(cover1diff ~ salt.timeseries, order = 2, data = bio1_salt)
grangertest(salt.timeseries ~ cover1diff, order = 2, data = bio1_salt)
grangertest(cover1diff ~ salt.timeseries, order = 3, data = bio1_salt)
grangertest(salt.timeseries ~ cover1diff, order = 3, data = bio1_salt)
grangertest(cover1diff ~ salt.timeseries, order = 4, data = bio1_salt)
grangertest(salt.timeseries ~ cover1diff, order = 4, data = bio1_salt)
grangertest(cover1diff ~ salt.timeseries, order = 5, data = bio1_salt)
grangertest(salt.timeseries ~ cover1diff, order = 5, data = bio1_salt)
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:1) + Lags(salt.timeseries, 1:1)
## Model 2: cover1diff ~ Lags(cover1diff, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.5277 0.4697
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:1) + Lags(cover1diff, 1:1)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:1)
##   Res.Df Df      F Pr(>F)
## 1     79                 
## 2     80 -1 0.5716 0.4519
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:2) + Lags(salt.timeseries, 1:2)
## Model 2: cover1diff ~ Lags(cover1diff, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.7189 0.4906
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:2) + Lags(cover1diff, 1:2)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:2)
##   Res.Df Df      F Pr(>F)
## 1     76                 
## 2     78 -2 0.8811 0.4185
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:3) + Lags(salt.timeseries, 1:3)
## Model 2: cover1diff ~ Lags(cover1diff, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.5766 0.2024
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:3) + Lags(cover1diff, 1:3)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:3)
##   Res.Df Df      F Pr(>F)
## 1     73                 
## 2     76 -3 1.2911 0.2839
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:4) + Lags(salt.timeseries, 1:4)
## Model 2: cover1diff ~ Lags(cover1diff, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 1.0591 0.3833
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:4) + Lags(cover1diff, 1:4)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:4)
##   Res.Df Df      F Pr(>F)
## 1     70                 
## 2     74 -4 0.9958 0.4157
## Granger causality test
## 
## Model 1: cover1diff ~ Lags(cover1diff, 1:5) + Lags(salt.timeseries, 1:5)
## Model 2: cover1diff ~ Lags(cover1diff, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 0.8786 0.5004
## Granger causality test
## 
## Model 1: salt.timeseries ~ Lags(salt.timeseries, 1:5) + Lags(cover1diff, 1:5)
## Model 2: salt.timeseries ~ Lags(salt.timeseries, 1:5)
##   Res.Df Df      F Pr(>F)
## 1     67                 
## 2     72 -5 0.6236 0.6822

No significant result.