Nachdem die Regelungen von Basel II jetzt einige Jahre im Einsatz sind, stellt sich mit der zunehmenden Datenhistorie die Frage, wie die Qualität der Prognosen gegen die Langzeitdurchschnitte zu validieren ist. Zum Thema der „Langzeithistorie der Jahresausfallraten“ mache ich hier mal einen Diskussionsvorschlag.
Am Ende des Artikels habe ich das Datenbeispiel und den Code beigefügt.
1. Prognose
Um das Beispiel besser nachvollziehbar zu machen, fange ich hier mit einer konstanten Prognose von jährlich 5,00% an:
PD Prognose 5% für alle Jahre
2. Jährlicher Binomialtest bei Irrtumswahrscheinlichkeit 10%
Aufgrund der unterschiedlichen Portfoliogröße in den Jahren sind die Konfidenzbänder unterschiedlich groß (rot: Ausfallrate wurde unterschätzt, grün: Ausfallrate wurde überschätzt).
Fehlerbalken Binomialtest 10% Irrtumswahrscheinlichkeit (beidseitig)
3. Vergleich mit Realisierung in jedem Jahr
Der Schätzung und den Konfidenzbändern aus dem Binomialtest werden dann die Realisierungen (hellgraue ausgefüllte Punkte) gegenübergestellt. Es ist zu sehen dass die Realisierung mal unter und mal über der Schätzung liegt – immer jedoch im Rahmen der jährlichen Fehlerspanne.
Jährliche Ausfallraten im Vergleich zur Prognose
4. Ermittlung der durchschnittlichen Jahresausfallraten
Der beobachtete Durchschnitt der Jahresausfallraten liegt relativ eng bei bei der jährlichen Prognose, jedoch im Laufe der 10 vorliegenden Jahre immer über der Prognose. Es stellt sich die Frage, ob eine signifikante Unterzeichnung des Risikos vorliegt.
Es wäre hilfreich, analog zum auf Jahresbasis durchgeführten Binomialtest ein Konfidenzband zu ermitteln, um daraus eine Teststatistik zu ermitteln.
Durchschnitt der Jahresausfallraten im Zeitablauf
5. Ermittlung der Verteilung aller Langzeithistorien – gegeben die PD-Prognose wäre zutreffend
Über die jährliche PD-Prognose und die Besetzung des Portfolios lässt sich über die Binomialverteilung mit Hilfe von R jeder denkbare Verlauf der Jahresausfallraten simulieren. Im Ergebnis zeigt sich die folgende Verteilung, aus der sich Konfidenzbänder ableiten lassen. Die Basis für den Test ist damit gelegt.
Simulierte Verteilung der durchschnittlichen Jahresausfallraten nach 10 Jahren
Auf jede der vorliegenden Perioden angewendet, ergibt sich ein durchaus plausibles Bild: Im Zeitablauf werden die aufgrund der Datenhistorie Konfidenzbänder enger. Je kürzer die Datenhistorie also ist, desto größere Abweichungen von der Prognose wären also zulässig.
Will man die Hypothese „Kalibrierung ist vor dem Hintergrund der Langfristhistorie korrekt.“ bei einem Alpha von 10% testen, so kann im speziellen Fall die Hypothese nicht widerlegt werden.
Konfidenzbänder um den Durchschnitt der Jahresausfallraten
Noch interessanter wird das Ganze, wenn die PD-Prognose von Jahr zu Jahr abweicht, in der Simulation ist das aber auch kein Problem.
6. Zusammenführung der Analysen
Die kombinierten Analysen auf 1-Jahres-Ebene mit der Langfristhistorie sind im letzten Chart aufgeführt:
Zusammenfassung aller Analysen
Fazit:
Trotz der Abweichungen im einjährigen Bereich und der leichten Unterschätzungen im bisher vorliegenden Langzeitdurchschnitt kann die Hypothese der korrekten Kalibriertheit bei einem Alpha von 10% nicht abgelehnt werden.
So einen Satz muss man einfach mal so stehen lassen, denke ich…
7. Code und Zahlenbeispiel:
Daten: Beispieldatei Langfristhistorie
Und hier kommt noch der Code:
###########################################################
#
#Idee zur Konstruktion eines statistischen Tests fuer den
#Langzeitdurchschnitt der Jahresausfallraten
#
###########################################################
# Generelle Beschreibung
# 1. Einlesen der vorhanden Historie (Schaetzung und
# Realisierung je Verfahren.
# 2. Simulieren der moeglichen Realisierungen, gegeben
# die Annahme der Binomialverteilung sowie PD-Prognose
# und Anzahl der Schuldner je Jahr
# 3. Ermittlung der Wahrscheinlichkeisverteilung fuer jede
# Realisierung der Langfristhistorie.
# 4. Test der Hypothese "Verfahren ist passend kalibriert"
# unter Beachtung der vorgegeben Irrtumswahrscheinlichkeiten
# (10%)
###########################################################
###########################################################
#Parametrisierung
###########################################################
#Anzahl Simulationen definieren
my.sim<-100000
#Irrtumswahrscheinlichkeiten definieren
my.alpha.rot<-0.10
#Pfad angeben
my.pfad<-paste("/Users/Christian/Documents/R-Projekte/",
"Langzeithistorie Ausfallraten/", sep="")
###########################################################
#Vorbereitung
###########################################################
#Arbeitsverzeichnis setzen
setwd(my.pfad)
#Historie einlesen
#Spalten: Jahr, Verfahren, Kunden, PD und Ausfaelle
#Wichtig:
# * Alle Verfahren muessen gleich viele Jahre haben.
# * Die Datensaetze je Verfahren folgen einander.
my.zeitreihe<-read.csv2("Langzeitdurchschnitt_Beispiel.csv",
header = TRUE, sep=";")
attach(my.zeitreihe)
#Vorhandene Verfahren bestimmen
my.verfahren<-levels(Verfahren)
#Anzahl Jahre ermitteln
my.jahre<-length(Jahr[Verfahren==my.verfahren[1]])
###########################################################
#Simulation
###########################################################
#Simulieren
#Ergebnis ist Matrix aus einer Spalte je Zeile der Historie und
# einer Zeile je Simulation
my.sim.ausfallraten function(i){rbinom(n=my.sim,
size=Kunden[i],
prob=PD[i])/Kunden[i]})
#Durchschnitt der Jahresausfallraten fuer jedes Verfahren fuer
#jede Simulation bestimmen
#Ergebnis ist eine Matrix aus einer Spalte je Verfahren und
# einer Zeile je Simulation
my.index<-1:my.jahre
my.durchschnitte.sim<-sapply(1:length(my.verfahren),
function(x){apply(my.sim.ausfallraten[,my.index+
(my.jahre*(x-1))], 1, mean)})
#Realisierte Durchschnitte berechnen
#Ergebnis ist ein Vektor mit einem Eintrag je Verfahren
my.durchschnitte.real<-sapply(1:length(my.verfahren),
function(x){mean(Ausfaelle[my.index+(my.jahre*(x-1))]/
Kunden[my.index+(my.jahre*(x-1))])})
#Schwellenwerte bestimmen
#Ergebnis ist jeweils ein Vektor mit einem Eintrag je Verfahren
my.schwellenwerte.gelb<-sapply(1:length(my.verfahren),
function(x){quantile(my.durchschnitte.sim[,x],
1-my.alpha.gelb)})
my.schwellenwerte.rot<-sapply(1:length(my.verfahren),
function(x){quantile(my.durchschnitte.sim[,x],
1-my.alpha.rot)})
#Test durchfuehren
#Ergebnis ist ein Vektor mit einem Eintrag je Verfahren
my.ampelfarbe<-array( ifelse(my.durchschnitte.real>my.schwellenwerte.rot,"rot",
ifelse(my.durchschnitte.real>my.schwellenwerte.gelb,
"gelb","gruen")))
###########################################################
#Ergebnisse aufbereiten
###########################################################
#Verteilung ausgeben
#Ergebnis ist eine Matrix aus einer Zeile fuer jedes Percentil
# und einer Spalte fuer jedes Verfahren
my.percentile<-c(seq(from=0, to=0.9, by=0.1),0.99, 1)
my.verteilung<-sapply(1:length(my.verfahren),
function(x){t(t(quantile(my.durchschnitte.sim[,x]
,my.percentile)))})
#Ergebnisse in Liste aufbereiten
my.erg function(x){list(name=Verfahren[1+(my.jahre*(x-1))]
, ampel=c(my.ampelfarbe[x])
, realisiert=c(my.durchschnitte.real[x])
, verteilung=matrix(data=my.verteilung[,x]
, nrow=length(my.verteilung[,x])
, ncol=1
, dimnames=list(perc=my.percentile,
c("Werte"))))})
#Grafiken erzeugen
my.plot<-function(x) {
my.max<-max(my.durchschnitte.sim[,x]
, my.durchschnitte.sim[,x]
, my.erg[,x]$realisiert)*1.05
my.min<-min(my.durchschnitte.sim[,x]
, my.durchschnitte.sim[,x]
, my.erg[,x]$realisiert)*0.95
hist(my.durchschnitte.sim[,x],
nclass=40,
main=paste("Simulierte Verteilung Durchschnitte der",
"Jahresausfallraten",
Verfahren[1+(my.jahre*(x-1))]),
xlim=c(my.min, my.max),
xlab="Durchnschnitt der Jahresausfallraten",
cex.main=0.8,
cex.axis=0.6,
cex.lab=0.5)
abline(v=my.erg[,x]$realisiert,col="grey",
lty="solid",lwd=3) #Realisierung
#abline(v=my.schwellenwerte.gelb[x],col="orange",lty=1,lwd=3)
abline(v=my.schwellenwerte.rot[x],col="red",lty="dotted",lwd=3)
}
for (i in 1:length(my.verfahren)) {
my.plot(i)
}
###########################################################
# Aufbereitung Zeitreihe
###########################################################
my.plot.zeitreihe<-function(x) {
#Prognose:
my.prog.mean<-PD[Verfahren==Verfahren[1+(my.jahre*(x-1))]]
my.prog.up<-qbinom(p=1-my.alpha.rot,
size=Kunden[Verfahren==Verfahren[1+(my.jahre*(x-1))]],
prob=PD[Verfahren==Verfahren[1+(my.jahre*(x-1))]])/
Kunden[Verfahren==Verfahren[1+(my.jahre*(x-1))]]
my.prog.down<-qbinom(p=my.alpha.rot,
size=Kunden[Verfahren==Verfahren[1+(my.jahre*(x-1))]],
prob=PD[Verfahren==Verfahren[1+(my.jahre*(x-1))]])/
Kunden[Verfahren==Verfahren[1+(my.jahre*(x-1))]]
my.prog.mean.avg<-sapply(1:length(my.prog.mean),
function(x){mean(my.prog.mean[1:x])})
my.index<-1:my.jahre+(my.jahre*(x-1))
my.durchschnitte.sim.lfd<-cbind(my.sim.ausfallraten[,my.index[1]],
sapply(my.index[2:length(my.index)],
function(i){apply(my.sim.ausfallraten[,my.index[1]:i],
1,mean)}))
my.prog.mean.up<-sapply(1:my.jahre,
function(i){quantile(my.durchschnitte.sim.lfd[,i],
1-my.alpha.rot)})
my.prog.mean.down<-sapply(1:my.jahre,
function(i){quantile(my.durchschnitte.sim.lfd[,i],
my.alpha.rot)})
#Realisierung
my.ausfallrate.1y<-Ausfaelle[Verfahren==Verfahren[1+(my.jahre*(x-1))]]/
Kunden[Verfahren==Verfahren[1+(my.jahre*(x-1))]]
my.real.mean.avg<-sapply(1:length(my.ausfallrate.1y),
function(x){mean(my.ausfallrate.1y[1:x])})
my.max<-max(my.prog.down,
my.prog.mean,
my.prog.mean.avg,
my.prog.mean.up,
my.prog.mean.down,
my.ausfallrate.1y,
my.real.mean.avg)*1.05
my.min<-min(my.prog.down,
my.prog.mean,
my.prog.mean.avg,
my.prog.mean.up,
my.prog.mean.down,
my.ausfallrate.1y,
my.real.mean.avg)*0.95
plot(my.prog.mean, type="p", pch=21 , col="grey", bg=NULL
, ylim=c(my.min, my.max)
, main=paste("Jahresprognose, Langzeithistorie\nund -realisierung",
Verfahren[1+(my.jahre*(x-1))])
, xlab="Jahr"
, ylab="Ausfallprognose / -rate"
, cex.main=0.8
, cex.axis=0.6
, cex.lab=0.5
)
arrows(y0=my.prog.mean, x0=1:length(my.prog.mean)
, y1=my.prog.up, length= 0.05, angle=90, col="red")
arrows(y0=my.prog.mean, x0=1:length(my.prog.mean)
, y1=my.prog.down, length= 0.05, angle=90, col="green")
lines(my.prog.mean.avg, lty="dotted", col="grey", lwd=2.5)
lines(my.prog.mean.up, lty="dotted", col="red", lwd=2.5)
lines(my.prog.mean.down, lty="dotted", col="green", lwd=2.5)
points(my.ausfallrate.1y, pch=21, col="grey", bg="grey", lwd=2.5)
lines(my.real.mean.avg, lty="solid", col="grey", lwd=2.5)
legend(x="topright", cex=0.5, bty="n", y.intersp=0.75,
c("Prog. mittl. Durchschn. 95%",
"Prog. mittl. Durchschn. 50%",
"Prog. mittl. Durchschn. 5%",
"Real. mittl. Durchschn.",
"Prog. 1-Jahres-PD",
"Real. 1-Jahres-Ausfallrate"),
lty=c("dotted", "dotted", "dotted", "solid", NA, NA, NA, NA),
lwd=c(2.5),
pch=c(NA, NA, NA, NA,21,21),
col=c("red", "grey", "green", "grey", "grey", "grey"),
pt.bg=c(rep(NA,5),"grey")
)
}
for (i in 1:length(my.verfahren)) {
#plot.new()
my.plot.zeitreihe(i)
}