Verinin yüklenmesi:

veri <-read.csv("https://raw.githubusercontent.com/gungorrbaris/zaman-serisi-analizi-R/main/data/Table_7.6_Electricity_End_Use.csv")
knitr::kable(head(veri,n=10), align = "c")

Month	Electricity.End.Use..Total
1973 January	144505.2
1973 February	139546.1
1973 March	137102.3
1973 April	131365.9
1973 May	131360.9
1973 June	140293.0
1973 July	152562.2
1973 August	157480.8
1973 September	157309.6
1973 October	146042.0

Verinin Düzenlenmesi:

Tarih değişkenini Ay-Yıl formatına getirelim:

tarih = seq(from = as.Date("1973-01-01"), to = as.Date("2021-11-01"), by = 'month')
veri$Month <- tarih
veri$Month <-format(veri$Month, "%m-%Y")
knitr::kable(head(veri), align = "c")

Month	Electricity.End.Use..Total
01-1973	144505.2
02-1973	139546.1
03-1973	137102.3
04-1973	131365.9
05-1973	131360.9
06-1973	140293.0

Değişken isimlerini türkçeleştirelim:

colnames(veri) <- c("Tarih","Elektrik_Tuketimi")
knitr::kable(head(veri), align = "c")

Tarih	Elektrik_Tuketimi
01-1973	144505.2
02-1973	139546.1
03-1973	137102.3
04-1973	131365.9
05-1973	131360.9
06-1973	140293.0

Zaman Serisinin Oluşturulması

veri_ts <- ts(data = veri$Elektrik_Tuketimi, start = c(1973, 01), end = c(2021, 11), frequency = 12 )

Verimizde gözlemler aylık olduğu için frekans değeri (frequency) 12 olarak belirlenmelidir.

Zaman Serisi Grafiğinin Çizilmesi

ts.plot(veri_ts,gpars=list(xlab="Zaman", ylab="Elektrik Tüketimi"))

Zaman Serisi Grafiğine bakıldığında:

• Pozitif yönlü artış vardır. Bu yüzden Serinin artan bir trende sahip olduğunu söyleyebiliriz.

• Seride düzenli dalgalanmaların olduğu görülmektedir. Bu seride mevsimsellik olabileceğini gösterir. Daha kesin sonuç için ACF grafiğini çizelim.

ACF Grafiğinin Çizilmesi

! BİLGİ !

Gecikmeli seri: Serilerdeki verilerin dönem kaydırılması işlemi ile ortaya çıkan serilerdir. Zt bir zaman serisi değişkeni ise; bir dönem gecikmeli zaman serisi Zt-1, iki dönem gecikmeli Zt-2 ve k dönem gecikmeli Zt-k ile gösterilir.

Otokorelasyon katsayısı: Zaman serileri ile bu serilerin gecikmeli serileri arasındaki ilişkileri verir.

ACF: Otokorelasyon fonksiyonu (Autocorrelation function) demektir. Otokorelasyon katsayısı değerlerinden oluşur.

Seride trendin olduğunu anlamak için ACF grafiğindeki ilk dört gecikmeye bakmak yeterlidir. İlk dört gecikme sınırlar dışındaysa seri için “trende sahiptir” diyebiliriz.

library(fpp)

## Zorunlu paket yükleniyor: forecast

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

## Zorunlu paket yükleniyor: fma

## Zorunlu paket yükleniyor: expsmooth

## Zorunlu paket yükleniyor: lmtest

## Zorunlu paket yükleniyor: zoo

## 
## Attaching package: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

## Zorunlu paket yükleniyor: tseries

Acf(veri_ts,lag.max = 42,  ylim=c(-1,1), lwd=5,main="Orijinal veri ACF Grafiği")

ACF Grafiğine bakıldığında:

• Tüm gecikmeler (4 gecikme olması yeterlidir) sınırlar dışında olduğu için zaman serisinin trende sahip olduğu söylenebilir. Bu yüzden serinin farkı alınıp durağan hale getirilmelidir.

• Trend ile birlikte düzenli dalgalanmaların olduğu gözükse de baskın bir mevsimsellikten şu an için bahsedemeyiz.

Birinci Farkın Alınması

! BİLGİ !

Fark Alınması: Zaman serisinin akışkanlı bir şekilde son değerlerinden belli bir dönem önceki değerlerinin çıkarılması işlemidir. Bu işlem sayesinde serideki trend ya da mevsimsel dalgalanmaları yok etmek mümkün olmaktadır.

veri_birinci_fark <- diff(veri_ts, differences = 1)

Farkı alınmış seri için ACF:

Acf(veri_birinci_fark,lag.max = 42,  ylim=c(-1,1), lwd=5,main="Birinci Fark sonrası ACF Grafiği")

MEVSİMSELLİK

Birinci Farkı alınmış seri için ACF Grafiğine bakıldığında:

• Yine tüm gecikmeler sınırlar dışında olduğu için birinci farkı alınmış zaman serisinin de trende sahip olduğu söylenebilir. Fakat bir öncekine göre bu ACF grafiğinde baskın bir mevsimsellik görülüyor. Dalgalanmalar ve gecikmelerdeki düzenli sıçramalar bunu destekliyor.

REGRESYON ANALİZİ

Zaman Serileri için Regresyon Analizinde, serinin mevsimsel olup olmamasına göre yöntemler değişiklik göstermektedir.

</blockquote>

Sadece trende sahip seriler (Mevsimsel olmayan) için;
Basit Doğrusal Regresyon
Birinci Farklar Regresyon Modeli
Üstel Regresyon Modeli
Karesel Regresyon Modeli
Lojistik Regresyon Modeli
Kübik Regresyon Modeli
gibi yöntemler kullanılabilir.

Hem trende hem de mevsimselliğe sahip seriler için ise önceki yazıda gördüğümüz Çarpımsal Ayrıştırma Yönteminden yararlanarak Regresyon analizini uygulayabiliriz.

Kaynak: Prof. Dr. Cem Kadılar, Dr. Hatice Öncel Çekim SPSS ve R Uygulamalı Zaman Serileri Analizine Giriş

Buradaki başlıkta; Trende sahip serinin birinci farkını aldığımızda serinin aynı zamanda mevsimselliğe sahip olduğunu ACF grafiğine bakarak açıklamıştık.

Kullandığımız veri mevsimselliğe sahip olduğu için bu bölümde sadece Mevsimsel Serilerde Regresyon Analizi yöntemleri gösterilecektir.

Birinci Harmonik

Oluşturulan yeni veri ile regresyon analizi yapılır.

regresyon.modelC1<-lm(data = veri_mevC_yeni, y ~ t+s1+c1)
summary(regresyon.modelC1)

## 
## Call:
## lm(formula = y ~ t + s1 + c1, data = veri_mevC_yeni)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -73643 -16474  -2108  16541  67571 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 149621.728   2066.679  72.397  < 2e-16 ***
## t              367.845      6.091  60.394  < 2e-16 ***
## s1             -46.962      4.296 -10.931  < 2e-16 ***
## c1             -27.767      4.306  -6.448 2.38e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 25000 on 583 degrees of freedom
## Multiple R-squared:  0.8679, Adjusted R-squared:  0.8672 
## F-statistic:  1276 on 3 and 583 DF,  p-value: < 2.2e-16

Birinci Harmonik için Regresyon modeli sonucuna göre:

• Coefficients kısmındaki p değerlerinin hepsi < α =0.05 olduğundan, modeldeki tüm terimlerin önemli olduğunu söyleyebiliriz.

• p-value değeri < α =0.05 olduğundan kurulan regresyon modelinin anlamlı olduğunu söyleyebiliriz.

İkinci Harmonik

Birinci harmonikte tüm terimler anlamlı çıktığından ikinci harmonikte de regresyon denklemine katılırlar.

Sin, Cos Terimlerinin oluşturulması:

s2<-t*sin(2*3.1416*2*t/12)
c2<-t*cos(2*3.1416*2*t/12)

Yeni Verinin Oluşturulması:

Yeni veri oluşturulurken birinci harmonikte kullanılan sin ve cos değerleri de eklenir.

veri_mevC_yeni2 <- as.data.frame(cbind(veri_ts, t, s1, c1, s2, c2))

names(veri_mevC_yeni2)<- c("y", "t", "s1", "c1", "s2", "c2")

knitr::kable(head(veri_mevC_yeni2), align = "c")

y	t	s1	c1	s2	c2
144505.2	1	0.5000011	0.8660248	0.8660266	0.4999979
139546.1	2	1.7320533	0.9999958	1.7320459	-1.0000085
137102.3	3	3.0000000	-0.0000110	-0.0000220	-3.0000000
131365.9	4	3.4640918	-2.0000170	-3.4641212	-1.9999661
131360.9	5	2.4999735	-4.3301423	-4.3300964	2.5000530
140293.0	6	-0.0000441	-6.0000000	0.0000882	6.0000000

İkinci Harmonik için regresyon modelinin kurulması:

regresyon.modelC2<-lm(data= veri_mevC_yeni2, y ~ t+s1+c1+s2+c2)
summary(regresyon.modelC2)

## 
## Call:
## lm(formula = y ~ t + s1 + c1 + s2 + c2, data = veri_mevC_yeni2)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -58872 -13619   -946  14783  45690 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 149492.533   1536.901  97.269  < 2e-16 ***
## t              368.503      4.529  81.356  < 2e-16 ***
## s1             -47.011      3.195 -14.715  < 2e-16 ***
## c1             -26.885      3.203  -8.394 3.58e-16 ***
## s2              69.020      3.195  21.605  < 2e-16 ***
## c2               8.234      3.203   2.571   0.0104 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 18590 on 581 degrees of freedom
## Multiple R-squared:  0.9272, Adjusted R-squared:  0.9266 
## F-statistic:  1479 on 5 and 581 DF,  p-value: < 2.2e-16

İkinci Harmonik için Regresyon modeli sonucuna göre:

• Coefficients kısmındaki p değerlerinin hepsi < α =0.05 olduğundan, modeldeki tüm terimlerin önemli olduğunu söyleyebiliriz.

• p-value değeri < α =0.05 olduğundan kurulan regresyon modelinin anlamlı olduğunu söyleyebiliriz.

Üçüncü Harmonik

İkinci harmonikte tüm terimler anlamlı çıktığından ikinci harmonikte de regresyon denklemine katılırlar.

Sin, Cos Terimlerinin oluşturulması:

s3<-t*sin(2*3.1416*3*t/12)
c3<-t*cos(2*3.1416*3*t/12)

Yeni Verinin Oluşturulması:

Yeni veri oluşturulurken birinci harmonikte kullanılan sin ve cos değerleri de eklenir.

veri_mevC_yeni3 <- as.data.frame(cbind(veri_ts, t, s1, c1, s2, c2,s3,c3))

names(veri_mevC_yeni3)<- c("y", "t", "s1", "c1", "s2", "c2","s3", "c3")

knitr::kable(head(veri_mevC_yeni3), align = "c")

y	t	s1	c1	s2	c2	s3	c3
144505.2	1	0.5000011	0.8660248	0.8660266	0.4999979	1.0000000	-3.70e-06
139546.1	2	1.7320533	0.9999958	1.7320459	-1.0000085	-0.0000147	-2.00e+00
137102.3	3	3.0000000	-0.0000110	-0.0000220	-3.0000000	-3.0000000	3.31e-05
131365.9	4	3.4640918	-2.0000170	-3.4641212	-1.9999661	0.0000588	4.00e+00
131360.9	5	2.4999735	-4.3301423	-4.3300964	2.5000530	5.0000000	-9.18e-05
140293.0	6	-0.0000441	-6.0000000	0.0000882	6.0000000	-0.0001322	-6.00e+00

Üçüncü Harmonik için regresyon modelinin kurulması:

regresyon.modelC3<-lm(data= veri_mevC_yeni3, y ~ t+s1+c1+s2+c2+s3+c3)
summary(regresyon.modelC3)

## 
## Call:
## lm(formula = y ~ t + s1 + c1 + s2 + c2 + s3 + c3, data = veri_mevC_yeni3)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -52489 -13350   -588  14875  44894 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 149478.431   1520.254  98.325  < 2e-16 ***
## t              368.575      4.480  82.263  < 2e-16 ***
## s1             -47.045      3.160 -14.887  < 2e-16 ***
## c1             -26.811      3.168  -8.463  < 2e-16 ***
## s2              68.916      3.160  21.808  < 2e-16 ***
## c2               8.319      3.168   2.626 0.008869 ** 
## s3               2.350      3.160   0.744 0.457328    
## c3              11.957      3.168   3.774 0.000177 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 18390 on 579 degrees of freedom
## Multiple R-squared:  0.929,  Adjusted R-squared:  0.9281 
## F-statistic:  1082 on 7 and 579 DF,  p-value: < 2.2e-16

Üçüncü Harmonik için Regresyon modeli sonucuna göre:

• Coefficients kısmındaki s3 terimi dışındaki p değerleri hepsi < α =0.05 olduğundan, s3 dışındaki terimlerin önemli olduğunu söyleyebiliriz.

• p-value değeri < α =0.05 olduğundan kurulan regresyon modelinin anlamlı olduğunu söyleyebiliriz.

• s3 terimi katsayısı = 0.457328 > α =0.05 olduğu için, s3 terimine ait katsayının anlamsız olduğunu söyleyebiliriz.

Bu durumda s3 değişkeni çıkarılarak tekrardan model anlamlılığına bakılır.

s3 teriminin olmadığı modelin oluşturulması:

regresyon.modelC4<-lm(data= veri_mevC_yeni3, y ~ t+s1+c1+s2+c2+c3)
summary(regresyon.modelC4)

## 
## Call:
## lm(formula = y ~ t + s1 + c1 + s2 + c2 + c3, data = veri_mevC_yeni3)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -52490 -13334   -518  14688  45695 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 149480.840   1519.665  98.364  < 2e-16 ***
## t              368.563      4.479  82.292  < 2e-16 ***
## s1             -47.045      3.159 -14.893  < 2e-16 ***
## c1             -26.825      3.167  -8.471  < 2e-16 ***
## s2              68.915      3.159  21.816  < 2e-16 ***
## c2               8.295      3.167   2.619 0.009038 ** 
## c3              11.957      3.167   3.776 0.000176 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 18380 on 580 degrees of freedom
## Multiple R-squared:  0.9289, Adjusted R-squared:  0.9282 
## F-statistic:  1263 on 6 and 580 DF,  p-value: < 2.2e-16

Modelin Anlamlılığı:

Ho: Model Anlamsızdır.
H1: Model Anlamlıdır.

• p-value değeri < α =0.05 olduğundan Ho RED. Kurulan regresyon modelinin anlamlı olduğunu söyleyebiliriz.

Durbin-Watson Testi:

Ho: Otokorelasyon sorunu yoktur.
H1: Otokorelasyon sorunu vardır.

library(fpp)

dwtest(data= veri_mevC_yeni3, y ~ t+s1+c1+s2+c2+c3)

## 
##  Durbin-Watson test
## 
## data:  y ~ t + s1 + c1 + s2 + c2 + c3
## DW = 0.2914, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

Durbin-Watson Test sonucuna göre:

• DW= 0,2914 değeri 2 ye yakın bir değer olmadığı için (p = 0 < α =0.05) Ho RED. Otokorelasyon sorunu vardır.

Çarpımsal Model ile kurulan Regresyon Analizinin Denklemi :

Zt = 149480.840 + 368.563t - 47.045sin(2πt/12) - 26.825cos(2πt/12) + 68.915sin(4πt/12) + 8.295cos(4πt/12) + 11.957cos(6πt/12) + εi olur.

Tahmin:

tahmin_modelC4<-predict(regresyon.modelC4)
plot.ts(tahmin_modelC4)

Alt ve Üst Sınırlar:

sinir_modelC4<-predict(regresyon.modelC4, interval = 'confidence' ,level = .95)
plot.ts(sinir_modelC4)

Hata:

hata_modelC4<-resid(regresyon.modelC4)
plot.ts(hata_modelC4)

Gerçek Elektrik Tüketimi ve Tahmin Sınırlarının Birlikte Çizilmesi

plot( window(veri_mevC_yeni3$y),
      xlab="", ylab="", type="l", lty=3, col=2, lwd=2)
lines(window(sinir_modelC4[,2]) ,type="l",lty=1,col=4,lwd=2)
lines(window(sinir_modelC4[,3]) ,type="l",lty=1,col=3,lwd=2)
legend("topleft",c(expression(paste(Elektri_Tüketimi)),
                   expression(paste(Altsinir)),
                   expression(paste(Üstsinir))),
       lwd=c(2,2,3),lty=c(3,1,2), cex=0.7, col=c(2,4,3))

plot( window(veri_mevC_yeni3$y), 
      xlab="Zaman", ylab="",type="l",lty=1, col=4, lwd=2)
lines( window(tahmin_modelC4),lty=3,col=2,lwd=3)
legend("topleft",c(expression(paste(Elektrik_tüketimi)),
                   expression(paste(Tahmin))),
       lwd=c(2,2),lty=c(1,3), cex=0.6, col=c(4,2))

Hatalar Akgürültü mü?

! BİLGİ !

Akgürültü Serisi: Durağanlık koşullarından tek farkı kovaryans teriminin sıfır olmasıdır. Dolayısıyla, akgürültü serisi durağan bir seriden farklı özellikler gösterir.
Örneğin, akgürültü serisi rasgele hareketlere sahip modellenemez bir seri iken durağan serilerin hareketlerinin belli bir sistematiği vardır ve bu nedenle modellenebilmektedir. Akgürültü serisinin tüm gecikmelerindeki otokorelasyon ve kısmi otokorelasyon değerleri önemsizdir.

ACF ve PACF grafıklerinin yorumunda serinin akgürültü serisi olup olmadığına net bir şekilde karar verilemiyorsa, yani ACF veya PACF grafiklerinde birinci gecikme dışında az sayıda güven sınırını biraz geçen ilişkiler var ise bu durumda seriye Box-Ljung Testi uygulanır.

Eğer her gecikme için Box-Ljung Testi sonucunda Ho : rk = 0 yokluk hipotezi kabul edilirse serinin akgürültü serisi olduğu söylenir.

Kaynak: Prof. Dr. Cem Kadılar, Dr. Hatice Öncel Çekim SPSS ve R Uygulamalı Zaman Serileri Analizine Giriş

Ho: Hatalar Akgürültüdür.(Hatalar arasında ilişki yoktur.)
H1: Hatalar Akgürültü değildir.(Hatalar arasında ilişki vardır.)

Hatalar için ACF Grafiği:

Acf(hata_modelC4, lag.max = 42,  ylim=c(-1,1), lwd=5,main="Hatalar için ACF Grafiği")

Hatalar için ACF Grafiğine bakıldığında:

• Tüm gecikmelerin sınırları geçtiği gözükmektedir.

• Hataların akgürültü olmadığı (Ho RED) Söylenebilir. Box-Ljung testi ile de bunu görebiliriz.

Box-Ljung Testi:

Box.test(hata_modelC4,type = "Ljung",lag = 42)

## 
##  Box-Ljung test
## 
## data:  hata_modelC4
## X-squared = 11835, df = 42, p-value < 2.2e-16

Box-Ljung test sonucuna bakıldığında:

• p değeri = 0 < α =0.05 ’dir. Yani Ho RED.

• Hatalar Arasında ilişki olduğunu yani hataların akgürültü olmadığını söyleyebiliriz.