This document presents the solution to the practical actuarial modelling exam.
library(fGarch)
## Warning: package 'fGarch' was built under R version 4.4.3
## NOTE: Packages 'fBasics', 'timeDate', and 'timeSeries' are no longer
## attached to the search() path when 'fGarch' is attached.
##
## If needed attach them yourself in your R script by e.g.,
## require("timeSeries")
library(timeSeries)
## Warning: package 'timeSeries' was built under R version 4.4.3
## Loading required package: timeDate
## Warning: package 'timeDate' was built under R version 4.4.3
##
## Attaching package: 'timeSeries'
## The following objects are masked from 'package:graphics':
##
## lines, points
import= read.csv("DATA25.csv")
D2=import$STOCK
plot(D2)
ts.plot(D2)
log.D2=diff(D2)
plot(log.D2)
ts.plot(log.D2)
acf(log.D2)
pacf(log.D2)
s=log.D2^2
plot(s)
ts.plot(s)
M2=garchFit(~garch(1,1),data = log.D2,trace = F)
summary(M2)
##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~garch(1, 1), data = log.D2, trace = F)
##
## Mean and Variance Equation:
## data ~ garch(1, 1)
## <environment: 0x000001a1bb87eb88>
## [data = log.D2]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## mu omega alpha1 beta1
## -0.0224632 0.0096866 0.8708402 0.3139453
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu -0.022463 0.023956 -0.938 0.34840
## omega 0.009687 0.005584 1.735 0.08279 .
## alpha1 0.870840 0.321604 2.708 0.00677 **
## beta1 0.313945 0.116863 2.686 0.00722 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## -11.23655 normalized: -0.09522499
##
## Description:
## Wed May 28 07:58:00 2025 by user: MAX
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 52.2611358 4.483747e-12
## Shapiro-Wilk Test R W 0.9321416 1.534461e-05
## Ljung-Box Test R Q(10) 23.4110865 9.326786e-03
## Ljung-Box Test R Q(15) 31.8778266 6.689206e-03
## Ljung-Box Test R Q(20) 38.5065827 7.674429e-03
## Ljung-Box Test R^2 Q(10) 3.2880888 9.738088e-01
## Ljung-Box Test R^2 Q(15) 4.5347064 9.953867e-01
## Ljung-Box Test R^2 Q(20) 6.7619136 9.973959e-01
## LM Arch Test R TR^2 6.4689585 8.906240e-01
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## 0.2582466 0.3521681 0.2560473 0.2963815
plot(M2, which=1)
plot(M2, which=6)
plot(M2,which=13)
##EGARCH Model
library(rugarch)
## Warning: package 'rugarch' was built under R version 4.4.3
## Loading required package: parallel
##
## Attaching package: 'rugarch'
## The following object is masked from 'package:stats':
##
## sigma
library(forecast)
## Warning: package 'forecast' was built under R version 4.4.3
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
C1= read.csv("DATA25.csv")
C2=C1$STOCK
C3=diff(D2) # Note: this uses D2 which was defined in q1_analysis chunk
M3=ugarchspec(variance.model=list(model="eGARCH",garchOrder=c(1,1)),mean.model=list(armaOrder=c(0,0),include.mean=F),distribution.model="norm")
M3fit=ugarchfit(C3,spec=M3)
plot(M3fit, which=10)
plot(M3fit, which=8)
plot(M3fit, which=1)
summary(M3fit, which=10)
## Length Class Mode
## 1 uGARCHfit S4
summary(M3fit,which=8)
## Length Class Mode
## 1 uGARCHfit S4
summary(M3fit,which=1)
## Length Class Mode
## 1 uGARCHfit S4
#plotting the forecast
forc=ugarchforecast(fitORspec=M3fit,n.ahead = 12)
plot(forc, which=1)
plot(forc, which=3)
##QUESTION TWO: Kaplan Meier models
library(survival)
library(survminer)
## Loading required package: ggplot2
## Loading required package: ggpubr
##
## Attaching package: 'ggpubr'
## The following object is masked from 'package:forecast':
##
## gghistogram
##
## Attaching package: 'survminer'
## The following object is masked from 'package:survival':
##
## myeloma
d1=data.frame(time=c(2,4,6,10),event=c(1,0,1,1))
kmc=with(d1,Surv(time,event))
plot(kmc)
plot(kmc,fun="cumhaz")
kmc2=survfit(Surv(time,event)~1,data=d1)
summary(kmc2)
## Call: survfit(formula = Surv(time, event) ~ 1, data = d1)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 2 4 1 0.750 0.217 0.4259 1
## 6 2 1 0.375 0.286 0.0839 1
## 10 1 1 0.000 NaN NA NA
##Kaplan-Meier with different data
library(survival)
library(survminer) # already loaded above, but harmless to include again
D1=data.frame(time=c(1,2,2,3,4,5,6,8,9,10),event=c(0,1,0,1,0,1,0,1,0,1))
kmc=with(D1,Surv(time,event))
kmc2 = surv_fit(Surv(time, event) ~ 1, data = D1)
summary(kmc2)
## Call: survfit(formula = Surv(time, event) ~ 1, data = structure(list(
## time = c(1, 2, 2, 3, 4, 5, 6, 8, 9, 10), event = c(0, 1,
## 0, 1, 0, 1, 0, 1, 0, 1)), class = "data.frame", row.names = c(NA,
## -10L)))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 2 9 1 0.889 0.105 0.706 1
## 3 7 1 0.762 0.148 0.521 1
## 5 5 1 0.610 0.181 0.341 1
## 8 3 1 0.406 0.205 0.151 1
## 10 1 1 0.000 NaN NA NA
plot(kmc,
col = "blue",
lwd = 2,
xlab = "Time",
ylab = "Survival Probability",
main = "Kaplan-Meier Curve")
grid() # Grid AFTER the plot call
plot(
kmc2,
fun = "cumhaz",
col = "purple",
lwd = 2,
xlab = "Time (Days)",
ylab = "Cumulative Hazard",
main = "Cumulative Hazard Function"
)
grid()