#--------------------------#
#                          #  
#   Loesungen Ueblatt 6    #
#                          #
#--------------------------#



# Wir benutzen die Funktionen
#
#  dDayMean()
#  dDayStdDev()
#  dDayNormRet()
#  
# aus dem week6.txt: 


dDayMean = function( d , ret )
{
  n = length(ret)
  result = rep(0,n)
  summe = 0

  for(i in 1:n)
  {
    summe = summe + ret[i]
    if(i > d) 
    {
      summe = summe - ret[i-d]
      result[i] = summe/d
    }
    else
    {
      result[i] = summe/i
    }
  }
  return(result)
}

dDayStdDev = function( d , ret )
{
  n = length(ret)
  result = rep(0,n)
  summe = 0

  for(i in 1:n)
  {
    summe = summe + ret[i]*ret[i]
    if(i > d) 
    {
      summe = summe - ret[i-d]*ret[i-d]
      result[i] = summe/d
    }
    else
    {
      result[i] = summe/i
    }
  }
  result = sqrt(result)    # square root element by element
  return(result)
}

dDayNormRet = function( d , ret )
{
  stddev = dDayStdDev(d,ret)

  # we have to shift by one element:
  stddev = stddev[-length(stddev)]   
  stddev = c(sd(ret),stddev)

  # stddev could be 0 if data are constant:
  stddev = pmax(stddev,0.00000001)     # "parallel max": maximum element by element

  result = ret/stddev                  # division element by element

  # ignore first d elements:
  result = result[-(1:d)]
  return(result)
}



#-------------------#
#  Analysis S&P500  #
#-------------------#

# Jetzt koennen wir die Daten laden und analysieren: 



spx = read.table("C:/Users/detlef/OneDrive/Hochschule/Vorlesungen/SS2025/Datenanalyse-mit-R/SPX.txt",header=TRUE,sep=";")
head(spx)
tail(spx)
summary(spx)
str(spx)                      # Adj.Close is numeric, no Factor

names(spx)
S = spx[,2]                   # ist aequivalent zu S = spx$Adj.Close
plot(S,type="l")
plot(log(S/S[1]),type="l")    # straight line would be constant exponential growth
n=length(S)
n

# wir berechnen die returns:

ret = rep(0,n)       
for(i in 2:n)
{
  ret[i] = (S[i]-S[i-1])/S[i-1]
}
plot(ret,type="l")


# let's take a look at Stylized Fact 1:

mean20 = dDayMean(20,ret)
mean60 = dDayMean(60,ret)
mean180 = dDayMean(180,ret)
plot(mean20)
points(mean60,col="red")
points(mean180,col="yellow")


# stddev's = daily volatilities:

stddev20 = dDayStdDev(20,ret)
stddev60 = dDayStdDev(60,ret)
stddev180 = dDayStdDev(180,ret)
stddev_n = dDayStdDev(n,ret)
plot(stddev20)
points(stddev60,col="red")
points(stddev180,col="yellow")
points(stddev_n,col="green")


# normalized returns:
normret15 = dDayNormRet(15,ret)
normret60 = dDayNormRet(60,ret)
normret250 = dDayNormRet(250,ret)
normret_n = ret/stddev_n[n]

par(mfrow=c(2,2))            # set up 2 time 2 plot array
                             # that is, 4 pictures at once

hist(normret_n,breaks=150,xlim=c(-5,5),ylim=c(0,0.8),prob=TRUE)
curve(dnorm(x),from=-5,to=5,add=TRUE,col="red")
hist(normret250,breaks=150,xlim=c(-5,5),ylim=c(0,0.8),prob=TRUE)
curve(dnorm(x),from=-5,to=5,add=TRUE,col="red")
hist(normret60,breaks=120,xlim=c(-5,5),ylim=c(0,0.8),prob=TRUE)
curve(dnorm(x),from=-5,to=5,add=TRUE,col="red")
hist(normret15,breaks=100,xlim=c(-5,5),ylim=c(0,0.8),prob=TRUE)
curve(dnorm(x),from=-5,to=5,add=TRUE,col="red")

# thus: if returns are normalized with more recent volatility 
# data, then normalized returns are more Gaussian 



#-----------------------------#
#  Analysis General Electric  #
#-----------------------------#

GE = read.table("C:/Users/detlef/OneDrive/Hochschule/Vorlesungen/SS2025/Datenanalyse-mit-R/GE.txt",header=TRUE,sep=";")
head(GE)
tail(GE)
summary(GE)                 # no NA's
names(GE)
S = GE[,2]                  
plot(S,type="l")
plot(log(S/S[1]),type="l")  # straight line would be constant exponential growth
n=length(S)
n
# jetzt berechnen wir die returns: -> redo above S&P500 code, conclusion also holds.