A Plex for all Seasons?
Aoife Fraoch
Monday, October 12, 2015
I would like to thank https://forums.eveonline.com/profile/Teckos%20Pech for the data used in this. The following is an exploration of PLEX pricing trends and seasonality. I am going to leave a reasonable portion of the code in this file visible in the interests of showing my working.
library(ggplot2)
library(reshape2)
library(forecast)## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## Loading required package: timeDate
## This is forecast 6.1library(data.table)##
## Attaching package: 'data.table'
##
## The following objects are masked from 'package:reshape2':
##
## dcast, meltlibrary(lmtest)
workingDirectory <- "C://data//other-projects//data" # Place all data files here #
dataFile <- "evePlexPrice.csv"
## Setting up the working directory for all data file references.
# Ensure all data files for use are in the working directory.
setwd(workingDirectory)
work_dir <- getwd()
## Function for exporting tables to CSVs.
##Directory path other than work directory and extension needs to be defined 'x', 'y' is the table to be exported.
file_output <- function (x, y){
path <- paste(work_dir, x, sep ="")
write.table(y, file = path, sep = ",", row.names = FALSE)
}
## Load data
evePlex.workFile <- read.csv(paste(workingDirectory,"//", dataFile, sep = ""), as.is=TRUE, sep="\t", skip = 1, quote='"', fill=TRUE, flush=TRUE)
## Clean up data.
evePlex.workFile[,"Average"] <- gsub("ISK","", evePlex.workFile[,"Average"])
evePlex.workFile[,"Average"] <- gsub(",","", evePlex.workFile[,"Average"])
evePlex.workFile[,"Average"] <- as.numeric(evePlex.workFile[,"Average"])
evePlex.workFile[,"Date"] <- as.Date(evePlex.workFile[,"Date"], "%m/%d/%Y")
evePlex.workFile[,"Quantity"] <- gsub(",","", evePlex.workFile[,"Average"])
evePlex.workFile[,"Quantity"] <- as.numeric(evePlex.workFile[,"Average"])Initial exploratory graphs
Average plex price ploted as per the imported data.
## First round trends.
ggplot(evePlex.workFile, aes(x = Date, y = Average)) + geom_line()
Time Series Decomposition
The following is a series of decomositions of the data by different periods.
## Month and Year Time Series.
evePlex.workFile <- evePlex.workFile[order(as.Date(evePlex.workFile$Date, format="%Y-%m-%d")),]
evePlex.workFile.dt <- data.table(evePlex.workFile)
PlexAverage <- evePlex.workFile.dt[,"Average", with = FALSE]
PlexAverage <- as.vector(PlexAverage$Average)
## So lets start with just a year period.
PlexAverage.ts = ts(PlexAverage, start = c(2010, 01, 01), frequency=365.25)
## period decomposition.
PlexAverage.stl = stl(PlexAverage.ts, s.window = "periodic")
plot(PlexAverage.stl)
summary(PlexAverage.stl)## Call:
## stl(x = PlexAverage.ts, s.window = "periodic")
##
## Time.series components:
## seasonal trend remainder
## Min. :-35681354 Min. : 282434178 Min. :-95443988
## 1st Qu.:-28619742 1st Qu.: 373887082 1st Qu.: -3094708
## Median :-21791600 Median : 528008716 Median : 14974313
## Mean :-19737704 Mean : 552250321 Mean : 16935094
## 3rd Qu.: -9368566 3rd Qu.: 696355266 3rd Qu.: 34690270
## Max. : -1828822 Max. :1018411264 Max. :303071245
## IQR:
## STL.seasonal STL.trend STL.remainder data
## 19251176 322468184 37784978 339937499
## % 5.7 94.9 11.1 100.0
##
## Weights: all == 1
##
## Other components: List of 5
## $ win : Named num [1:3] 20901 549 365
## $ deg : Named int [1:3] 0 1 1
## $ jump : Named num [1:3] 2091 55 37
## $ inner: int 2
## $ outer: int 0The Durbin-Watson test does seem to indicate that there is more than just a 365.25 day seasonal period.
dwtest(tslm(PlexAverage.ts ~ trend + season), alt="two.sided")##
## Durbin-Watson test
##
## data: tslm(PlexAverage.ts ~ trend + season)
## DW = 0.028495, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is not 0Initial decomposition. While a seasonal and trend componant were both extracted, the residuals account for more variance than the seasonal componant.
## 7 period decomposition.
PlexAverageW.stl = stl(PlexAverage.ts, s.window = 7)
plot(PlexAverageW.stl)
summary(PlexAverageW.stl)## Call:
## stl(x = PlexAverage.ts, s.window = 7)
##
## Time.series components:
## seasonal trend remainder
## Min. :-41130929 Min. :274872658 Min. :-80360744
## 1st Qu.:-14398698 1st Qu.:377943594 1st Qu.:-18161745
## Median : -3504559 Median :526587111 Median : -2946289
## Mean : -383325 Mean :551684905 Mean : -1853869
## 3rd Qu.: 9265216 3rd Qu.:693054568 3rd Qu.: 13314787
## Max. : 92190956 Max. :992941481 Max. :203328766
## IQR:
## STL.seasonal STL.trend STL.remainder data
## 23663914 315110974 31476533 339937499
## % 7.0 92.7 9.3 100.0
##
## Weights: all == 1
##
## Other components: List of 5
## $ win : Named num [1:3] 7 699 365
## $ deg : Named int [1:3] 0 1 1
## $ jump : Named num [1:3] 1 70 37
## $ inner: int 2
## $ outer: int 0## 365.25 period decomposition.
PlexAverageY.stl = stl(PlexAverage.ts, s.window = 365.25)
plot(PlexAverageY.stl)
summary(PlexAverageY.stl)## Call:
## stl(x = PlexAverage.ts, s.window = 365.25)
##
## Time.series components:
## seasonal trend remainder
## Min. :-35530692 Min. : 282322716 Min. :-83921979
## 1st Qu.:-14724880 1st Qu.: 373883993 1st Qu.:-19159635
## Median : -6046214 Median : 528030195 Median : -2788057
## Mean : -1099763 Mean : 552272536 Mean : -1725063
## 3rd Qu.: 10259703 3rd Qu.: 696430145 3rd Qu.: 14363568
## Max. : 59883154 Max. :1018049509 Max. :222120914
## IQR:
## STL.seasonal STL.trend STL.remainder data
## 24984583 322546152 33523203 339937499
## % 7.3 94.9 9.9 100.0
##
## Weights: all == 1
##
## Other components: List of 5
## $ win : Named num [1:3] 365 551 365
## $ deg : Named int [1:3] 0 1 1
## $ jump : Named num [1:3] 37 56 37
## $ inner: int 2
## $ outer: int 0For the 7 and 365.25 day periods a seasonal and trend componant were also extracted. However in both of those cases the residuals account for more variance than the seasonal componant.
The residuals do appear to be quiet significant in the last year when compared to previous years in the data set. Especially during the last spike in plex prices.
Further Approaches
PlexAverage.tbats <- tbats(PlexAverage.ts, seasonal.periods = c(7, 365.25), use.box.cox = TRUE, use.arma.errors = TRUE)
plot(PlexAverage.tbats)
plot(PlexAverage.tbats$errors)
PlexAverage.tbats## TBATS(0, {0,0}, 0.818, {<7,2>, <365.25,6>})
##
## Call: tbats(y = PlexAverage.ts, use.box.cox = TRUE, seasonal.periods = c(7,
## 365.25), use.arma.errors = TRUE)
##
## Parameters
## Lambda: 1e-06
## Alpha: 0.8484456
## Beta: -0.0420649
## Damping Parameter: 0.818149
## Gamma-1 Values: -6.835356e-07 -0.0003498706
## Gamma-2 Values: -0.000415182 6.235535e-05
##
## Seed States:
## [,1]
## [1,] 19.3749637057
## [2,] 0.0011979326
## [3,] -0.0013182453
## [4,] 0.0001372439
## [5,] -0.0052594233
## [6,] -0.0020033266
## [7,] 0.0138799037
## [8,] -0.0051783423
## [9,] 0.0016202060
## [10,] -0.0040753685
## [11,] -0.0038637206
## [12,] -0.0058311206
## [13,] -0.0300494801
## [14,] -0.0257461946
## [15,] -0.0009936129
## [16,] 0.0025961112
## [17,] 0.0060952254
## [18,] 0.0071881899
##
## Sigma: 0.01311054
## AIC: 81755.66Additional work
From the looks of the data, there are a number of extreme shortlived events that should be removed from the model. Possibly a more stable model could be created if the data was aggregated from days.
There is definitely an upwards trend, and there also does appear to be a yearly trend. It seems that random activity is more influencial on day to day prices however.
Forcast Models
Create an stl forcast model based on the initial timeseries object.
PlexAverage.stlf <- stlf(PlexAverage.ts, lambda = TRUE, h = 90)
plot(PlexAverage.stl$time.series)
PlexAverage.stlf$model## ETS(M,N,N)
##
## Call:
## ets(y = x, model = etsmodel, allow.multiplicative.trend = allow.multiplicative.trend)
##
## Smoothing parameters:
## alpha = 0.8583
##
## Initial states:
## l = 267301205.1582
##
## sigma: 0.0118
##
## AIC AICc BIC
## 81275.77 81275.77 81287.06plot(PlexAverage.stlf$seasonal)
stl model
plot(PlexAverage.stlf)
tbats model
plot(forecast(PlexAverage.tbats, h = 90))