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

У зависности од тога како се подаци региструју временске серије могу бити прекидне и непрекидне. У зависности од тога да ли се статистичка својства серије мењају током времена серије могу бити стационарне и нестационарне.

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) # odredjuje klasu objekta 
## [1] "ts"
start(AP) # pocetak vs 
## [1] 1949    1
end(AP) # kraj vs 
## [1] 1960   12
frequency(AP) # period vs
## [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)")

\(\color{deepskyblue}{\text{Компоненте временске серије}}\)

plot(aggregate(AP)) # uklanjanje sezonske varijacije, izravnavanje serije, tako da jasnije vidimo trend

boxplot(AP~cycle(AP)) #sve vrednosti iz AP dele se u 12 kategorija i za svaku od njih se crta boxplot

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

Сезонска компонента AP - највише путника има лети, а најмање крајем године.


\(\color{deepskyblue}{\text{Циљеви анализе временске серије}}\)

  1. Описивање временске серије
  2. Избор одговарајућег математичког модела
  3. Предвиђање временске серије

\(\color{deepskyblue}{\text{Основни појмови у анализи временских серија}}\)


Maine <- read.csv("Maine.txt", sep="")
attach(Maine)
class(Maine) 
## [1] "data.frame"
Maine.month.ts <- ts(unemploy, start=c(1996,1), freq=12) # konvertovanje u ts
Maine.annual.ts <- aggregate(Maine.month.ts)/12 # pravimo sravnjenu seriju - od mesecne dobijamo godisnju
par(mfrow=c(1,2))
plot(Maine.month.ts)
plot(Maine.annual.ts)

(Maine.Feb <- window(Maine.month.ts, start = c(1996,2), freq = TRUE))
## Time Series:
## Start = 1996.083 
## End = 2006.083 
## Frequency = 1 
##  [1] 6.7 6.5 5.7 5.0 4.4 4.2 4.9 5.8 5.6 5.8 5.6
(Maine.Aug <- window(Maine.month.ts, start = c(1996,8), freq = TRUE))
## Time Series:
## Start = 1996.583 
## End = 2006.583 
## Frequency = 1 
##  [1] 4.0 4.0 3.6 3.3 2.5 3.1 3.6 4.3 3.8 4.1 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% мања од просека.


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

Видимо сличности између серија - опада стопа незапослености почетком ’00. Разлика је што у Maine-у стопа незапосленоти не опада 2006. као у остатку САД.


CBE <- read.delim("cbe.txt") # ovo je multiple vremenska serija
CBE[1:4, ] # prve 4 vrste i sve kolone
##   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)
Beer.ts <- ts(CBE[, 2], start = 1958, freq = 12)
Choc.ts <- ts(CBE[, 1], start = 1958, freq = 12)
par(mfrow=c(1,1))
plot(cbind(Elec.ts, Beer.ts, Choc.ts))# sa funkcijom cbind unutar plot, dobijamo vise vremenskih serija

Примећујемо растући тренд у све три временске серије.

AP.elec <- ts.intersect(AP, Elec.ts) # dve serije spajamo u jednu dvodimenzionu na preseku vremena
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")

par(mfrow=c(1,1))
plot(as.vector(AP),as.vector(Elec), xlab ="Air passengers/1000's", ylab="Electricity production / MWh") 
#potrebno je pretvoriti serije u nizove
abline(reg = lm(Elec ~ AP)) 

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

Видимо да су линеарно зависне, преко корелације која мери линеарну зависност. Међутим, иако су корелисане, не значи да су једна последица друге (AP у USA, Elec у AUS)!


Z <- pounds_nz <- read.csv("pounds_nz.txt", sep="")# exchange rate - koliko jedne valute se dobija za jedinicnu vrednost druge (r), ili obrnuto (bice 1/r)
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)) #nismo naveli frequency=TRUE, pa uzima sve vrednosti u navedenom intervalu
par(mfrow=c(1,2))
plot(Z.92.96)
plot(Z.96.98)

par(mfrow=c(1,1))

global <- scan("global.txt",  sep="") # ucita u niz
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) # najvise nas zanima trend pa izravnjujemo seriju
plot(Global.ts)

plot(Global.annual)

New.series <- window(Global.ts, start=c(1970,1)) #izdvajamo samo onaj period kritican za globalno zagrevanje
plot(New.series)
New.time <- time(New.series) # funkcija time u niz stavlja sve vremenske trenutke u kojima su vrsena merenja
abline(lm(New.series~New.time)) # malopre smo imali regresiju kojom smo ispitivali odnos dva obelezja,а sada trend ocenjujemo regresionom pravom

\(\color{deepskyblue}{\text{Општи модели временских серија}}\)

Адитивни модел: \(x_t= m_t + s_t + z_t,\)

Мултипликативни модел: \(x_t=m_t*s_t*z_t\)

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

\(\color{deepskyblue}{\text{Како оценити тренд?}}\)

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

Како изабрати \(a_i\)? Специјално, \(a_i=\frac{a}{2m+1}\).

Ако претпоставимо да је сезона једна година, да би се анулирала сезонска компонента треба сабрати свих 12 месеци. Међутим, ако оцењујемо тренд у јулу требало би да се одлучимо да ли ћемо сабрати 5 месеци пре јула и 6 после или обрнуто. Да би се избегла замена, узима се средина израчунатих вредности у оба случаја. Одавде се добије оцена сезонске компоненте одузимањем, односно дељењем почетне серије.

\(\color{deepskyblue}{\text{Декомпозиција временске серије}}\)

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=2)

\(\color{deepskyblue}{\text{Функције за оцењивање тренда и сезонске компоненте методом покретних просека}}\)

Направимо функцију за оцењивање тренда методом покретних просека. Резултат треба да буде исти као тренд добијен функцијом \(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(Trend)

plot(MA(Elec.ts))

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

За детрендовану временску серију се израчунају средње вредности по сваком месецу, па се тај низ од 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()\) за адитивни модел. Слично је и за мултипликативни модел.

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