# Opens and assign variables.

#setwd(choose.dir())

master.raw <-read.csv("master.csv")


head(master.raw)
##   year        pdo       enso   enso_pdo       fire   melt
## 1 1984  0.8375000 -0.6508333 -0.5450729  10.991895 121.85
## 2 1985  0.4491667 -0.6400000 -0.2874667  18.991054 113.95
## 3 1986  1.2391667  0.1566667  0.1941361   0.000000 108.90
## 4 1987  1.8208333  1.2958333  2.3594965   3.614698 111.20
## 5 1988  0.5316667 -0.9483333 -0.5041972 446.850532 109.55
## 6 1989 -0.1791667 -0.7183333  0.1287014   4.039639 109.45
fire.raw <- master.raw$fire # [1000 HA's]
snow.raw <- master.raw$melt # [Avg melt date since Jan 1st]
pdo.raw <- master.raw$pdo # [PDO index]
enso.raw <- master.raw$enso # [ENSO Index]
enso.pdo.raw <- master.raw$enso_pdo # [ ENSO * PDO]
years <- master.raw$year # [Years]

Normality Testing of Raw Data

# Test for normality using shapiro.test

shapiro.test(fire.raw) # Data is not normal
## 
##  Shapiro-Wilk normality test
## 
## data:  fire.raw
## W = 0.592, p-value = 8.239e-08
shapiro.test(snow.raw)
## 
##  Shapiro-Wilk normality test
## 
## data:  snow.raw
## W = 0.9732, p-value = 0.6487
shapiro.test(pdo.raw)
## 
##  Shapiro-Wilk normality test
## 
## data:  pdo.raw
## W = 0.9775, p-value = 0.7719
shapiro.test(enso.raw)
## 
##  Shapiro-Wilk normality test
## 
## data:  enso.raw
## W = 0.958, p-value = 0.2934
shapiro.test(enso.pdo.raw)
## 
##  Shapiro-Wilk normality test
## 
## data:  enso.pdo.raw
## W = 0.8365, p-value = 0.0004056

All data is normal except the fire data. One way to possibly make the data normal is to transform the data.

Transformations

Trying to make the fire data normal.

# Try to transform fire data to make it normal

# Square Root
fire.raw.sqrt <- sqrt(fire.raw+1)
shapiro.test(fire.raw.sqrt)
## 
##  Shapiro-Wilk normality test
## 
## data:  fire.raw.sqrt
## W = 0.7947, p-value = 6.719e-05
# Log
fire.raw.log <- log(fire.raw+100)
shapiro.test(fire.raw.log)
## 
##  Shapiro-Wilk normality test
## 
## data:  fire.raw.log
## W = 0.7184, p-value = 3.904e-06
# 1/x
fire.raw.recip <- 1/(fire.raw+1)
shapiro.test(fire.raw.recip)
## 
##  Shapiro-Wilk normality test
## 
## data:  fire.raw.recip
## W = 0.6222, p-value = 1.919e-07
# Cube Root
fire.raw.cuberoot <- fire.raw^1/3
shapiro.test(fire.raw.cuberoot)
## 
##  Shapiro-Wilk normality test
## 
## data:  fire.raw.cuberoot
## W = 0.592, p-value = 8.239e-08

None of the transformations have made these data normal. Thus we canot use pearson correlation, and instead must spearman's correlation to test for significance, since the fire data is not normal.

Z-scores and probability

# Convert to Z-scores using scale function. 

fire.z <- scale(fire.raw)
snow.z <- scale(snow.raw)
enso.z <- scale(enso.raw)
pdo.z <- scale(pdo.raw)
enso.pdo.z <- scale(enso.pdo.raw)

Using Z-scores should not yield a different p-value than the raw data, but I tried it anyway.

Testing correlations

# Evaluate correlations using three methods, Pearson being the preferred method for these data.
cor.test(fire.z, snow.z, method="kendall")
## Warning in cor.test.default(fire.z, snow.z, method = "kendall"): Cannot
## compute exact p-value with ties
## 
##  Kendall's rank correlation tau
## 
## data:  fire.z and snow.z
## z = -1.8389, p-value = 0.06592
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## -0.2419753
cor.test(fire.z, snow.z, method="spearman")
## Warning in cor.test.default(fire.z, snow.z, method = "spearman"): Cannot
## compute exact p-value with ties
## 
##  Spearman's rank correlation rho
## 
## data:  fire.z and snow.z
## S = 5402.831, p-value = 0.07969
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.3307465
cor.test(fire.z, snow.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z and snow.z
## t = -1.2967, df = 27, p-value = 0.2057
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.5590256  0.1364892
## sample estimates:
##        cor 
## -0.2421287
## PDO/ENSO - Fire comparisons
cor.test(fire.z, pdo.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z and pdo.z
## t = -1.0056, df = 27, p-value = 0.3236
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.5202713  0.1897227
## sample estimates:
##        cor 
## -0.1899941
cor.test(fire.z, enso.pdo.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z and enso.pdo.z
## t = -1.2753, df = 27, p-value = 0.2131
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.5562620  0.1404231
## sample estimates:
##        cor 
## -0.2383494
cor.test(fire.z, enso.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z and enso.z
## t = -1.9814, df = 27, p-value = 0.05782
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.63931421  0.01174754
## sample estimates:
##        cor 
## -0.3562922
#  Plot the fire/snow quadrants.
cor.test(snow.z, pdo.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  snow.z and pdo.z
## t = 1.468, df = 27, p-value = 0.1537
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1051045  0.5805335
## sample estimates:
##       cor 
## 0.2718741
cor.test(snow.z, enso.pdo.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  snow.z and enso.pdo.z
## t = 1.0623, df = 27, p-value = 0.2975
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1793776  0.5280387
## sample estimates:
##       cor 
## 0.2002961
cor.test(snow.z, enso.z, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  snow.z and enso.z
## t = 0.3043, df = 27, p-value = 0.7632
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.3147906  0.4160511
## sample estimates:
##        cor 
## 0.05845969
plot(fire.z, snow.z, xlab="Z-value Total Area Burned", ylab="Z-value Snow")
abline(h=0)
abline(v=0)

Without 1988

It looks like it doesn’t really matter if you re-scale the data from the raw file, the correlation test results come out the same.

#firedata<-read.csv("S:/Huxley/ENVS-497E_Pyro_Medler/R_Statistics/CSVs/Fire_Snowmelt_Zscores_Final_no_1988.csv")

fire.raw88=fire.raw[c(-5)]
snow.raw88=snow.raw[c(-5)]
fire.z88=scale(fire.raw88)
snow.z88=scale(snow.raw88)

# Evaluate correlations using three methods, Pearson being the preferred method for these data.
cor.test(fire.z88, snow.z88, method="kendall")
## Warning in cor.test.default(fire.z88, snow.z88, method = "kendall"):
## Cannot compute exact p-value with ties
## 
##  Kendall's rank correlation tau
## 
## data:  fire.z88 and snow.z88
## z = -2.0555, p-value = 0.03983
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## -0.2758621
cor.test(fire.z88, snow.z88, method="spearman")
## Warning in cor.test.default(fire.z88, snow.z88, method = "spearman"):
## Cannot compute exact p-value with ties
## 
##  Spearman's rank correlation rho
## 
## data:  fire.z88 and snow.z88
## S = 5021.874, p-value = 0.04969
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.3743498
cor.test(fire.z88, snow.z88, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z88 and snow.z88
## t = -2.461, df = 26, p-value = 0.02081
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.69502828 -0.07349944
## sample estimates:
##        cor 
## -0.4346576
#  Plot the fire/snow quadrants.
plot(fire.z88, snow.z88, xlab="Z-value Total Area Burned", ylab="Z-value Snow")
abline(h=0)
abline(v=0)

Most Recent 20 years 1993-2012

fire.z93=fire.z[c(-1,-2,-3,-4,-5,-6,-7,-8,-9)]
snow.z93=snow.z[c(-1,-2,-3,-4,-5,-6,-7,-8,-9)]

# Evaluate correlations using three methods, Pearson being the preferred method for these data.
cor.test(fire.z93, snow.z93, method="kendall")
## Warning in cor.test.default(fire.z93, snow.z93, method = "kendall"):
## Cannot compute exact p-value with ties
## 
##  Kendall's rank correlation tau
## 
## data:  fire.z93 and snow.z93
## z = -2.1749, p-value = 0.02964
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## -0.3535632
cor.test(fire.z93, snow.z93, method="spearman")
## Warning in cor.test.default(fire.z93, snow.z93, method = "spearman"):
## Cannot compute exact p-value with ties
## 
##  Spearman's rank correlation rho
## 
## data:  fire.z93 and snow.z93
## S = 1925.224, p-value = 0.04786
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.4475367
cor.test(fire.z93, snow.z93, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z93 and snow.z93
## t = -2.1386, df = 18, p-value = 0.04643
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.744369811 -0.009485938
## sample estimates:
##        cor 
## -0.4501173
#  Plot the fire/snow quadrants.
plot(fire.z93, snow.z93, xlab="Z-value Total Area Burned", ylab="Z-value Snow")
abline(h=0)
abline(v=0)

Most Recent 12 years 2000-2012

fire.raw2000=fire.raw[c(-16:-1)]
snow.raw2000=snow.raw[c(-16:-1)]
fire.raw2000
##  [1] 253.42814  45.45305  33.49471 109.21431   1.72166  12.98232 144.02255
##  [8] 103.67440  40.54913  14.75448  12.86172  38.06173 282.67130
snow.raw2000
##  [1]  90.15  96.60 106.25 100.30  83.75  79.15  94.05  84.60 110.85 112.30
## [11] 103.05 120.90  90.90
fire.z2000=scale(fire.raw2000)
snow.z2000=scale(snow.raw2000)

# Evaluate correlations using three methods, Pearson being the preferred method for these data.
cor.test(fire.z2000, snow.z2000, method="kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  fire.z2000 and snow.z2000
## T = 34, p-value = 0.59
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## -0.1282051
cor.test(fire.z2000, snow.z2000, method="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  fire.z2000 and snow.z2000
## S = 402, p-value = 0.737
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.1043956
cor.test(fire.z2000, snow.z2000, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z2000 and snow.z2000
## t = -1.0022, df = 11, p-value = 0.3378
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.7247433  0.3113372
## sample estimates:
##        cor 
## -0.2892702
#  Plot the fire/snow quadrants.
plot(fire.z2000, snow.z2000, xlab="Z-value Total Area Burned", ylab="Z-value Snow")
abline(h=0)
abline(v=0)

Playing with Splines

library(zoo)
## Warning: package 'zoo' was built under R version 3.1.3
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(lattice)
## Warning: package 'lattice' was built under R version 3.1.3
master.zoo=read.zoo(master.raw, format = "%Y" )
master.zoo
##                   pdo        enso    enso_pdo       fire   melt
## 1984-05-29  0.8375000 -0.65083333 -0.54507292  10.991895 121.85
## 1985-05-29  0.4491667 -0.64000000 -0.28746667  18.991054 113.95
## 1986-05-29  1.2391667  0.15666667  0.19413611   0.000000 108.90
## 1987-05-29  1.8208333  1.29583333  2.35949653   3.614698 111.20
## 1988-05-29  0.5316667 -0.94833333 -0.50419722 446.850532 109.55
## 1989-05-29 -0.1791667 -0.71833333  0.12870139   4.039639 109.45
## 1990-05-29 -0.3558333  0.24750000 -0.08806875  15.453414  93.30
## 1991-05-29 -0.4191667  0.66750000 -0.27979375  14.467016 107.65
## 1992-05-29  0.9283333  0.69416667  0.64441806   2.509013  93.85
## 1993-05-29  1.4166667  0.50666667  0.71777778   0.000000 117.15
## 1994-05-29 -0.1516667  0.45083333 -0.06837639  24.263203 104.00
## 1995-05-29  0.6425000  0.04583333  0.02944792   2.010305 112.50
## 1996-05-29  0.6408333 -0.29833333 -0.19118194  17.137333 111.80
## 1997-05-29  1.4608333  1.26250000  1.84430208   3.561563 127.50
## 1998-05-29  0.2458333  0.14333333  0.03523611   1.704471 117.15
## 1999-05-29 -1.0633333 -1.02500000  1.08991667   9.920696  99.05
## 2000-05-29 -0.5900000 -0.80333333  0.47396667 253.428136  90.15
## 2001-05-29 -0.5625000 -0.26750000  0.15046875  45.453046  96.60
## 2002-05-29  0.2208333  0.65916667  0.14556597  33.494712 106.25
## 2003-05-29  0.9691667  0.30666667  0.29721111 109.214310 100.30
## 2004-05-29  0.3450000  0.43916667  0.15151250   1.721660  83.75
## 2005-05-29  0.3750000  0.10583333  0.03968750  12.982320  79.15
## 2006-05-29  0.1908333  0.16500000  0.03148750 144.022550  94.05
## 2007-05-29 -0.1958333 -0.48916667  0.09579514 103.674399  84.60
## 2008-05-29 -1.2925000 -0.64166667  0.82935417  40.549128 110.85
## 2009-05-29 -0.6125000  0.42250000 -0.25878125  14.754475 112.30
## 2010-05-29 -0.3125000 -0.36166667  0.11302083  12.861723 103.05
## 2011-05-29 -1.2308333 -0.77750000  0.95697292  38.061732 120.90
## 2012-05-29 -1.1000000  0.03750000 -0.04125000 282.671300  90.90
xyplot(master.zoo)

fire.spline <- smooth.spline(master.raw$fire,spar=.5)$y
plot(master.raw$fire, type="l", col="blue")
lines(fire.spline)

melt.spline <- smooth.spline(master.raw$melt,spar=.5)$y
plot(master.raw$melt, type="l", col="blue")
lines(melt.spline)

melt.z.spline=scale(melt.spline)
fire.z.spline=scale(fire.spline)
cor.test(fire.z.spline, melt.z.spline, method="kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  fire.z.spline and melt.z.spline
## T = 146, p-value = 0.03315
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## -0.2807882
cor.test(fire.z.spline, melt.z.spline, method="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  fire.z.spline and melt.z.spline
## S = 5944, p-value = 0.01197
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.4640394
cor.test(fire.z.spline, melt.z.spline, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.z.spline and melt.z.spline
## t = -1.9557, df = 27, p-value = 0.06092
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.6365724  0.0163701
## sample estimates:
##       cor 
## -0.352249

Splines on the 2000-2012 data

melt.spline2000 <- smooth.spline(snow.raw2000,spar=.5)$y


fire.spline2000 <- smooth.spline(fire.raw2000,spar=.5)$y
plot(fire.raw2000, type="l", col="blue")
lines(fire.spline2000)

melt.spline2000 <- smooth.spline(snow.raw2000,spar=.5)$y
plot(snow.raw2000, type="l", col="blue")
lines(melt.spline2000)

melt.z.spline=scale(melt.spline)
fire.z.spline=scale(fire.spline)
cor.test(fire.spline2000, melt.spline2000, method="kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  fire.spline2000 and melt.spline2000
## T = 41, p-value = 0.8577
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## 0.05128205
cor.test(fire.spline2000, melt.spline2000, method="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  fire.spline2000 and melt.spline2000
## S = 354, p-value = 0.935
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## 0.02747253
cor.test(fire.spline2000, melt.spline2000, method="pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  fire.spline2000 and melt.spline2000
## t = 0.3379, df = 11, p-value = 0.7418
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.4762267  0.6178359
## sample estimates:
##      cor 
## 0.101353

ACF()

snow.ts=ts(master.raw$melt, start=1984)
plot(snow.ts)

snow.ts.acf=acf(snow.ts)

fire.ts=ts(master.raw$fire, start=1984)
plot(fire.ts)

fire.ts.acf=acf(fire.ts)