## Tarvitset paketin fGarch (Packages->install package(s) -> "paikka" --> fGarch)
## sekä tämän tarvitsevat paketetit timeDate, timeSeries ja fBasics; ks. ?fGarch

## Luetaan mainittuun hakemistoon tallennettu aikasarja - tässä tapauksessa Apple

data=read.table(file.choose(), header=TRUE)  ## Tämän komennon pitäisi antaa sinun valita kansioistasi haluamasi tiedosto

## Vaihtoehtoinen tapa. Muuta tiedoston polku oikeaksi ja poista kommentointi jos haluat käyttää tätä menetelmää.
##
## data=read.table("C:/Opiskelu/Tilastotiede/Stationaariset 2015/Laskuharjoitukset/Apple.txt",header=TRUE)

y = ts(data$Y)  ## muodostetaan aikasarja

## Aikasarjan piirtäminen (ks. R-koodi_1)
## Aikasarjan muunnokset (myös differenssointi) ja keskistys tarvittaessa (ks. R-koodi_1)

## GARCH -mallin estimointi; ?garchFit
## 'include.mean = FALSE', jos  käytetään keskistettyjä havaintoja tai odotusarvo oletetaan nollaksi; 
## 'include.mean = TRUE', jos mallissa estimoitava vakiotermi
## Mahdollisten jakaumien joukossa mm. "norm" (normaalijakauma), "std" (t-jakauma) ja "QMLE" (normaalijakauma 
## robustein keskivirhein) sekä vinot versiot kahdesta ensin mainitusta, "snorm" ja "sstd";   ks. ?fGarch.
## Huomaa, että garch( , ):in 1. argumentti on ARCH-aste eli toisin päin kuin kurssilla.

## Tavanomaisen GARCH-mallin estimointi (tässä virheen jakaumaksi valittu normaalijakauma robustein keskivirhein)

fit <- garchFit(~garch(s,r),y,include.mean = TRUE, include.shape = NULL, cond.dist = c("QMLE"))
fit

summary(fit)   ## hieman yksityiskohtaisempi tulostus, jossa myös Ljung-Box- ja McLeod-Li -testi (Jarque-Bera -testi on testi, 
              ## joka testaa virheiden normaalisuutta käyttäen standardoitujen residuaalien kolmansia ja neljänsiä momentteja)
plot(fit)   ## diagnostiikkaa graafisesti; sisältää mm. standardoitujen residuaalien ja niiden neliöiden otosautokorrelaatiofunktiot ja qq-plotit

## Estimoidun ehdollisen varianssin tai hajonnan aikasarjakuvan piirtäminen.
## 'type = c("sigma")', jos valitaan  ehdollinen hajnonta,  'type = c("h")', jos valitaan  ehdollinen varianssi

vola <- volatility(fit, type = c("sigma"))
vola <- ts(vola)
plot(vola)

## ARMA-GARCH -mallin estimointi (GARCH'in 1. argumentti tässäkin ARCH-aste eli toisin päin kuin kurssilla); ?garchFit
## Vakion käsittely ja jakaumat kuten GARCH-tapauksessa (tässä virheen jakaumaksi valittu t-jakauma)

fit <- garchFit(~arma(p,q)+garch(s,r),y,include.mean = TRUE, include.shape = NULL, cond.dist = c("std"))
fit

summary(fit) ## kuten edellä GARCH-tapauksessa
plot(fit)  ## kuten edellä GARCH-tapauksessa

vola <- volatility(fit, type = c("sigma"))  ## kuten edellä GARCH-tapauksessa
vola <- ts(vola)
plot(vola)