Временске серије и примене у финансијама

Стефан Малбашић, Математички факултет, Универзитет у Београду

Обавезе:

  • писмени део испита вреди 50 поена (за полагање је довољно освојити 20 поена)
  • усмени део испита вреди 50 поена

Литература:

  • Paul S.P. Cowpertwait · Andrew V. Metcalfe “Introductory Time Series with R”
  • Jonathan D. Cryer, Kung-Sik ChanTime “Series Analysis With Applications in R”
  • Peter J. Brockwell, Richard A. Davis “Introduction to Time Series and Forecasting”
  • Rob J Hyndman, George Athanasopoulos “Forecasting: Principles and Practice”

Основни појмови

Временска серија је низ података (опсервација) са дефинисаним редоследом, најчешће хронолошким. Временска серија може да се схвати и као случајни процес \(\{X(t), t\in T\}\), са непрекидним или дискретним параметром, при чему засеци \(X(t)\) и \(X(s)\), \(t\neq s\) у општем случају нису некорелисани. Анализа временских серија се разликује од стандардне теорије случајних процеса (у којој је полазни појам прост случајан узорак) и базира се на испитивању природе и степена корелисаности. Важне претпоставке које омогућују извођење статистичких закључака су стационарност и ергодичност временске серије.

Постоје различити типови временских серија, а поделу можемо вршити по следећим критеријумима:

  • непрекидне и прекидне серије;
  • једнодимензионе и вишедимензионе серије;
  • стационарне и нестационарне серије.

Air passenger bookings

Број међународних путничких резервација (у хиљадама) месечно на авионској компанији (Pan Am) у Сједињеним Америчким Државама добијен је од Федералне авијационе управе за период 1949-1960. године (Браун, 1963). Компанија је користила податке да предвиди будућу потражњу пре него што наручи нове авионе и обучи ваздухопловно особље. Подаци су доступни као временска серија у R и илуструју неколико важних концепата који се појављују у истраживачкој анализи временских серија.

data(AirPassengers)
AP <- AirPassengers
AP
##      Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1949 112 118 132 129 121 135 148 148 136 119 104 118
## 1950 115 126 141 135 125 149 170 170 158 133 114 140
## 1951 145 150 178 163 172 178 199 199 184 162 146 166
## 1952 171 180 193 181 183 218 230 242 209 191 172 194
## 1953 196 196 236 235 229 243 264 272 237 211 180 201
## 1954 204 188 235 227 234 264 302 293 259 229 203 229
## 1955 242 233 267 269 270 315 364 347 312 274 237 278
## 1956 284 277 317 313 318 374 413 405 355 306 271 306
## 1957 315 301 356 348 355 422 465 467 404 347 305 336
## 1958 340 318 362 348 363 435 491 505 404 359 310 337
## 1959 360 342 406 396 420 472 548 559 463 407 362 405
## 1960 417 391 419 461 472 535 622 606 508 461 390 432
class(AP) # одређује класу објекта 
## [1] "ts"
start(AP) # почетак серије
## [1] 1949    1
end(AP) # крај серије
## [1] 1960   12
frequency(AP) # период серије
## [1] 12
summary(AP)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   104.0   180.0   265.5   280.3   360.5   622.0
plot(AP, ylab = "Passengers (1000s)")

Компоненте временске серије

  • Тренд је дугорочна правилност у кретању временске серије. Може бити детерминистички (описује се неком детерминистичком функцијом, нпр. глобално загревање) и стохастички (који се не може предвидети на основу историјских података, односно реализоване серије, нпр. финансије).

  • Сезонска компонента обухвата правилност у кретању временске серије која се појављује у току једне календарске године. Већи је степен корелације између опсервација које одговарају истим месецима у различитим година, него различитим месецима у истој години (температура у неком месту).

  • Цикличне варијације (циклуси) су правилности у кретању временске серије које се појављују са периодом различитим (обично већим) од једне године. Потребно је одредити период цикличне варијације (пролазак комете).

  • Случајна компонента представља промене у кретању временске серије које су случајног карактера (нпр. поплаве или суше у серији која се односи на количину падавина у неком месту). То је једина компонента која обавезно постоји у временској серији, јер би, у супротном, био у питању детерминистички низ.

За добијање јаснијег увида у тренд, сезонска компонента се може уклонити агрегирањем података на годишњем нивоу, што се може постићи коришћењем функције aggregate(која од месечне серије прави годишњу). Преглед вредности за сваку сезону може се видети коришћењем boxplot-a, при чему се cycle функција користи за издвајање сезона за сваку ставку података.

plot(aggregate(AP))

boxplot(AP ~ cycle(AP)) # све вредности из серије AP се деле у 12 категорија и за сваку од њих се црта boxplot

Разлог за растући тренд AP је напредовање авионске индустрије, повећање популације, бољи животни стандард након Другог светског рата.

Сезонска компонента AP - највише путника има лети (јун-септембар).

Циљеви анализе временске серије

Техника рада са (реализованим) временским серијама састоји се од следећих корака:

  1. Описивање - припрема временске серије за даљу анализу,

  2. Објашњење - избор одговарајућег математичког модела по којем се мењају вредности серије,

  3. Предвиђање (прогнозирање) - на основу вредности серије у неком периоду у прошлости, треба предвидети њене вредности у будућности.

Основни појмови у анализи временских серија

  • Случајни процес

  • Временска серија

  • Строга и слаба стационарност

  • Бели шум

  • Ергодичност

Unemployment: Maine

Стопа незапослености један је од главних макроекономских показатеља. На пример, она утичу на политике за регионални развој и обезбеђивање економске стабилности. Месечна стопа незапослености за америчку државу Мејн од јануара 1996. до августа 2006. године приказана дата је у бази Maine.

Maine <- read.csv("Maine.txt", sep="")
attach(Maine)
class(Maine)
## [1] "data.frame"
Maine.month.ts <- ts(Maine, start=c(1996,1), freq=12) # конвертовање у објекат класе "ts"
Maine.annual.ts <- aggregate(Maine.month.ts)/12 # агрегирање серије (од месечне добијамо годишњу)
plot(Maine.month.ts)

plot(Maine.annual.ts)

Постоји јасна месечна варијација. Из претходних графика се види да је опсервација за фебруар обично отприлике за 20% већа од годишњег просека, док је опсервација за август обично отприлике за 20% мања.

(Maine.Feb <- window(Maine.month.ts, start = c(1996,2), freq = TRUE))
## Time Series:
## Start = 1996.083 
## End = 2006.083 
## Frequency = 1 
##       unemploy
##  [1,]      6.7
##  [2,]      6.5
##  [3,]      5.7
##  [4,]      5.0
##  [5,]      4.4
##  [6,]      4.2
##  [7,]      4.9
##  [8,]      5.8
##  [9,]      5.6
## [10,]      5.8
## [11,]      5.6
(Maine.Aug <- window(Maine.month.ts, start = c(1996,8), freq = TRUE))
## Time Series:
## Start = 1996.583 
## End = 2006.583 
## Frequency = 1 
##       unemploy
##  [1,]      4.0
##  [2,]      4.0
##  [3,]      3.6
##  [4,]      3.3
##  [5,]      2.5
##  [6,]      3.1
##  [7,]      3.6
##  [8,]      4.3
##  [9,]      3.8
## [10,]      4.1
## [11,]      3.9
(Feb.ratio <- mean(Maine.Feb)/mean(Maine.month.ts))
## [1] 1.222529
(Aug.ratio <- mean(Maine.Aug)/mean(Maine.month.ts))
## [1] 0.8163732

Стопа незапослености је 22% виша у фебруару од просечне и 18% нижа у августу. Објашњење за ово је да Мејн привлачи туристе током лета, што креира више радних места. Такође, период пре Божића и преко новогодишњих празника обично има више стопе запослености него првих неколико месеци нове године.

Unemployment: US

Месечна стопа незапослености за целе Сједињене Америчке Државе од јануара 1996. до октобра 2006. године приказана је на графику испод.

US.month <- read.csv("USunemp.txt", sep="") 
attach(US.month)
US.month.ts <- ts(US.month, start=c(1996,1), freq = 12)
plot(US.month.ts, ylab = "unemployed (%)")

Смањење стопе незапослености почетком 2000. године је заједничко за Мејн и целе Сједињене Америчке Државе, али Мејн се не чини да дели смањење незапослености 2004. године као у остатку Сједињених Америчких Држава.

Вишедимензиона временска серија

Овде илуструјемо неколико важних идеја и концепата повезаних са вишедимензионим временским серијама. Месечна испорука електричне енергије (у милионима kWh), пива (у милионима литара) и чоколадних производа (у тонама) у Аустралији у периоду од јануара 1958. до децембра 1990. године доступни су од стране Аустралијског бироа за статистику (ABS). Ове три серије су смештене у један датотеку, која се може учитати на следећи начин:

CBE <- read.delim("cbe.txt") 
CBE[1:4, ] # прве 4 врсте и све колоне
##   choc beer elec
## 1 1451 96.3 1497
## 2 2037 84.4 1463
## 3 2477 91.2 1648
## 4 2785 81.9 1595
class(CBE) 
## [1] "data.frame"
Elec.ts <- ts(CBE[, 3], start = 1958, freq = 12) # не мора да се наведе други аргумент за start или end, тада се претпоставља да су једнаки 1 odnosno 12, редом
Beer.ts <- ts(CBE[, 2], start = 1958, freq = 12)
Choc.ts <- ts(CBE[, 1], start = 1958, freq = 12)
plot(cbind(Elec.ts, Beer.ts, Choc.ts)) # са функцијом cbind унутар plot, добијамо више временских серија

Графикони показују растуће трендове у све три временске серије, делом због пораста популације у Аустралији са око 10 милиона на око 18 милиона током истог периода. Али приметимо да је производња електричне енергије порасла за око 7 пута, а производња чоколаде за око 4 пута током овог периода, а да се популација није сасвим удвостручила.

AP.elec <- ts.intersect(AP, Elec.ts) # две серије спајамо у једну дводимензиону на пресеку времена
start(AP.elec)
## [1] 1958    1
end(AP.elec)
## [1] 1960   12
AP.elec[1:3,]
##       AP Elec.ts
## [1,] 340    1497
## [2,] 318    1463
## [3,] 362    1648
AP <- AP.elec[,1] 
Elec <- AP.elec[,2] 
par(mfrow=c(1,2))
plot(AP, main = "", ylab = "Air passengers / 1000's")
plot(Elec, main = "", ylab = "Electricity production / MkWh")

Графици броја међународних путника и производње електричне енергије у Аустралији за период 1958-1960 изгледају слично зато што обе серије имају растући тренд и сезонску компоненту. Међутим, ово не подразумева да постоји узрочно-последични однос између променљивих.

par(mfrow=c(1,1))
plot(as.vector(AP),as.vector(Elec), xlab ="Air passengers/1000's", ylab="Electricity production / MWh") 
# потребно је претворити серије у низове
abline(reg = lm(Elec ~ AP)) 

cor(AP,Elec) 
## [1] 0.8841668

Временске серије су високо корелисане са коефицијентом корелације 0.88. Међутим, иако су корелисане, не значи да су једна последица друге (AP у USA, Elec у AUS)! Разумно објашњење за корелацију је технолошки развој у обе земље током овог периода који доводи до растућег тренда.

Структурни лом: GBP to NZ dollar

У претходна два примера, трендови и сезонске компоненте били су јасни из графика. Поред тога, могу се изнети разумна објашњења за могуће узроке ових карактеристика. Са економским подацима, као што је кретање дивизних курсева, такве изражене шеме мање су вероватне да ће бити видљиве, и обично су потребни различити методи анализе. Финансијска серија понекад може показати драматичну промену која има јасан узрок, као што су рат или природна катастрофа. Промене из дана у дан тешко је објаснити јер су подлежни условима који су комплексни и немогући за изоловање, и често ће бити нереалистично претпоставити било какве детерминистичке компоненте у моделу таквих временских серија.

Кретање девизног курса за британску фунту у односу на новозеландски долар за период од јануара 1991. до марта 2000. године приказано је у наставку. Подаци су просечне вредности за кварталне периоде од три месеца, први квартал је од јануара до марта, а последњи квартал је од октобра до децембара.

Z <- pounds_nz <- read.csv("pounds_nz.txt", sep="")
Z.ts <- ts(Z, frequency = 4, start = 1991)
plot(Z.ts, xlab = "time (years)", ylab = "quarterly exchange rate in $NZ/pound")

Можемо издвојити краткорочне (short-term) трендове, пад од 1992, до 1996. године, па раст од 1996. до 1998. године. Смер кретања тренда се мења (опадајући, растући) у непредвидим тренуцима. За овакву временску серију кажемо да има структурни лом.

Z.92.96<-window(Z.ts, start = 1992, end = c(1996,1))
Z.96.98<-window(Z.ts, start = 1996, end = c(1998,1)) # нисмо навели frequency=TRUE, па узима све вредности у наведеном интервалу
par(mfrow = c(1,2))
plot(Z.92.96)
plot(Z.96.98)

par(mfrow = c(1,1))

Global temperature series

У студијама о климатским променама (на пример, видети Џонс и Моберг, 2003; Рејнер и други, 2003), следећа серија глобалних температура, изражена као одступање од месечних средњих вредности за период од 1856. до 2005. године, игра централну улогу. Највише нас занима тренд, зато се функција aggregate користи за уклањање сезонских ефеката унутар сваке године и производи годишњу серију просечних температура. Можемо избећи експлицитно дељење са 12 ако наведемо FUN = mean у функцији aggregate.

global <- scan("global.txt",  sep="") # учитавање података у низ
Global.ts <- ts(global, frequency = 12, start = c(1856,1), end = c(2005,12))
start(Global.ts)
## [1] 1856    1
end(Global.ts)
## [1] 2005   12
Global.annual <- aggregate(Global.ts, FUN=mean) 
plot(Global.ts)

plot(Global.annual)

New.series <- window(Global.ts, start = c(1970,1)) # издвајамо само онај период критичан за глобално загревање 
plot(New.series)
New.time <- time(New.series) # функција time у низ ставља све временске тренутке у којима су вршена мерења 
abline(lm(New.series~New.time)) # тренд серије оцењујемо регресионом правом

Општи модели временских серија

Адитивни модел: \(x_t= f_t + s_t + z_t\) (све компоненте серије се сабирају)

Мултипликативни модел: \(x_t=f_t\cdot s_t\cdot z_t\) (све компоненте серије се множе)

  • \(f_t\) вредност тренда серије у тренутку \(t\)
  • \(s_t\) вредност сезонске компоненте у тренутку \(t\)
  • \(z_t\) вредност случајне компоненте серије у тренутку \(t\) (процес грешке), у општем случају корелисане случајне величине са очекивањем 0

Мултипликативни модел се користи ако сезонска компонента расте како расте тренд. Своди се на адитивни када логаритмујемо (серија мора бити позитивна). Проблем настаје када користимо функцију \(exp()\) да бисмо се вратили на оригинални модел. Грешка нема нормалну расподелу.

Како оценити тренд?

  • До сада смо покушали да оценимо тренд регресијом.

Метод покретних просека је непараметарски метод за изравнање временске серије, који за оцену тренда у тачки \(t\) користи просек неколико елемената у околини те тачке, за свако покретно \(t\). Посматрајмо низ реализованих опсервација: \[x_1,x_2,\dots,x_{t-m-1},x_{t-m},\dots,x_t,\dots,x_{t+m},x_{t+m+1},\dots,x_n.\] Сваку опсервацију \(x_t\) замељујемо са опсервацијом \(y_t\), где је \[y_t=\sum_{i=-m}^{m}a_i\cdot x_{t+i}, \quad a_i>0, \forall i, \quad \sum_{i=-m}^{m}a_i=1.\] Вредност тренда у тренутку \(t\) се оцењује са \(y_t\): \[\hat{f}_t = y_t.\] На овај начин се могу добити оцене тренда у тренуцима \(t=m+1,m+2,\dots,n-m\). Добијамо мању серију (одбацили смо првих \(m\) и последњих \(m\) чланова), која је глаткија. Наведени поступак може се применити више пута. За веће \(m\) добија се глаткија серија, али се избаци више чланова.

Како изабрати пондере \(a_i\)? Често се користи најједноставнији случај методе покретних просека, који се добија за \[a_i=\frac{1}{2m+1}, \quad t = -m,\dots,m.\]

Ако претпоставимо да је сезона једна година, да би се анулирала сезонска компонента треба сабрати свих 12 месеци. Међутим, ако оцењујемо тренд у јулу требало би да се одлучимо да ли ћемо сабрати 5 месеци пре јула и 6 после или обрнуто. Да би се избегла замена, узима се средина израчунатих вредности у оба случаја. \[{y}_t=\frac{\frac{1}{2}x_{t-6}+x_{t-5}+x_{t-4}+x_{t-3}+x_{t-2}+x_{t-1}+x_{t}+x_{t+1}+x_{t+2}+x_{t+3}+x_{t+4}+x_{t+5}+\frac{1}{2}x_{t+6}}{12}\]

Декомпозиција временске серије

  • Користи метод покретних просека за оцену тренда
Elec.decom.add <- decompose(Elec.ts)
plot(Elec.decom.add)

Elec.decom.mult <- decompose(Elec.ts, type = "mult")
plot(Elec.decom.mult)

На основу првог графика видимо да је погоднији мултипликативни модел, јер сезонска компонента расте са порастом тренда.

Trend <- Elec.decom.add$trend #nedostaje prvih 6 i poslednjih 6 elemenata
Seasonal <- Elec.decom.add$seasonal
ts.plot(cbind(Trend,Trend+Seasonal),lty=1:2)

Функције за оцењивање тренда и сезонске компоненте методом покретних просека

Направимо функцију за оцењивање тренда методом покретних просека. Резултат треба да буде исти као тренд добијен функцијом \(decompose()\).

MA<-function(ts){
  x<-vector()
  ts1<-as.vector(ts)
  x[1:6]<-NA; x[(length(ts1)-6):(length(ts1))]<-NA
  for(i in 7:(length(ts1)-6)) 
  {
    x[i]<-(ts1[i-6]/2+ts1[i+6]/2+sum(ts1[(i-5):(i+5)]))/12
  }
  x.ts<-ts(x,start=c(start(ts)[1],start(ts)[2]),freq=12)
  return(x.ts)
}
plot(MA(Elec.ts))

plot(Trend)

Још је потребно видети како се тачно добија сезонска компонента.

За детрендовану временску серију се израчунају средње вредности по сваком месецу, па се тај низ од 12 вредности центрира.

SV<-function(ts,trend){
  m<-t(matrix(data = ts-trend, nrow = 12))
  seasonal<-colMeans(m, na.rm = TRUE)
  seasonal<-seasonal-mean(seasonal)
  ts.seasonal<-ts(rep(seasonal,end(ts)[1]-start(ts)[1]+1),start=start(ts)[1],freq=12)
  plot(ts.seasonal,ylab="Seasonal Variation")
  return(seasonal)
}
SV(Elec.ts,Trend)

##  [1] -498.75966 -639.83127 -150.13205 -402.49403  277.95649  525.28331
##  [7]  920.25987  654.66482   53.53071  -27.88205 -318.18544 -394.41070
plot(Seasonal)

Добили смо два иста графика и уверили се како тачно ради функција \(decompose()\) за адитивни модел. Слично је и за мултипликативни модел.

Ручно смо добили тренд и сезонску компоненту, а грешку добијамо када од оригиналне серије одузмемо тренд и сезонску компоненту.