data("lynx")   #loads the lynx dataset
str(lynx)      #shows structure of dataset lynx

plot(lynx,ylab="annual number of lynx",xlab="time in years")
plot(log(lynx),ylab="log annual number of lynx",xlab="time in years")
plot(log(lynx),type="b",ylab="log annual number of lynx",xlab="time in years")

qqnorm(lynx)    #QQ normal plot
qqnorm(log(lynx))     #QQ normal plot of log-transformed data

data("sunspots")  #loads the sunspots dataset
str(sunspots)
plot(sunspots,main="monthly number of sunspots")

qqnorm(sunspots)
qqnorm(log(0.5+sunspots))                                        #not successful
qqnorm(sqrt(sunspots))            #not successful

plot(stl(sunspots,s.window="periodic"))
plot(stl(sunspots,s.window=121))
decomp <- stl(sunspots,s.window="periodic")
decomp$time.series[1:13,]

acf(sunspots)
acf(decomp$time.series[,3])

par(mfrow=c(1,2))
acf(sunspots)
acf(decomp$time.series[,3])

par(mfrow=c(2,2))
plot(sunspots)
plot(decomp$time.series[,3])
acf(sunspots)
acf(decomp$time.series[,3])

par(mfrow=c(1,1))
qqnorm(decomp$time.series[,3])    #marginal Gaussianity is OK

#ARCH(1) model

#ARCH-Prozesse
f.arch.sim <- function(n=500,init=200,a)
{
  ## Purpose: Simulating the Arch(1) process
  ## -------------------------------------------------------------------------
  ## Arguments: n=length of the series, init=length of the initial segment
  ## that is discarded, a= Parameter of the model
  ## -------------------------------------------------------------------------
  ## Author: Hans-Ruedi Kuensch, Date:  5 Feb 97, 16:49
m <- n+init
x <- rnorm(m)
for (i in (2:m)) x[i] <- sqrt(1+a*x[i-1]^2)*x[i]
ts(x[(init+1):m],start=1)}

set.seed(22)
x <- f.arch.sim(n=1000,a=0.9)

plot(x)
acf(x)

plot(x)
acf(x)
acf(x^2)
acf(abs(x))
acf(sqrt(abs(x)))
acf(sin(x))
acf(sin(6*x))
acf(x^3)
acf(sign(x))


#--> uncorrelatedness but not independent


#AR(p) model simluation

x <- arima.sim(n=1000,model=list(ar=c(0.8,-0.7)))
plot(x)
acf(x)

x <- arima.sim(n=1000,model=list(ar=c(0.5,0.3)))
#does not work: non-stationary model



#Simulation of random-walk
help(filter)   #gives details about the function "filter"

set.seed(22)
x <- filter(ts(rnorm(1000)),filter=1,method="recursive",init=rnorm(1))
                                #generates X[t] = X[t-1] + eps[t]
plot(x)
acf(x)
plot(diff(x))
acf(diff(x))