---
title:
output:
  revealjs::revealjs_presentation:
    theme: white
    css: white-edits.css
    highlight: tango
    transition: slide
    self_contained: false
    reveal_plugins: ["chalkboard","search"]
    transitionSpead: fast
    mathjax: "https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.1/MathJax.js?config=TeX-AMS_SVG"
    reveal_options:
      minScale: 1.0
      maxScale: 1.0
      width: 1536
      height: 864
      slideNumber: true
      controlsTutorial: true
      chalkboard:
        theme: whiteboard
        readOnly: false
---

## ACCT 420: Logistic Regression for Bankruptcy {data-background-video="../Backgrounds/title.mp4" data-background-video-loop=true class="header no-toc-progress present-only" data-state="no-toc-progress"}

<!-- "reveal_plugins:..." and "chalkboard:" are handled by prep38.py -->
```{r helpers, echo=FALSE, warning=FALSE, message=F}
library(knitr)
library(kableExtra)
html_df <- function(text, cols=NULL, col1=FALSE, full=F) {
  if(!length(cols)) {
    cols=colnames(text)
  }
  if(!col1) {
    kable(text,"html", col.names = cols, align = c("l",rep('c',length(cols)-1))) %>%
      kable_styling(bootstrap_options = c("striped","hover"), full_width=full)
  } else {
    kable(text,"html", col.names = cols, align = c("l",rep('c',length(cols)-1))) %>%
      kable_styling(bootstrap_options = c("striped","hover"), full_width=full) %>%
      column_spec(1,bold=T)
  }
}
```

<h3 class="date"><br>Session 5</h3>
<h2 class="author"><br>Dr. Richard M. Crowley<br>rcrowley@smu.edu.sg<br><a href="http://rmc.link/">http://rmc.link/</a></h2>

## ACCT 420: Logistic Regression for Bankruptcy {data-background="../Backgrounds/title.jpg" class="header no-toc-progress present-not" data-state="no-toc-progress"}

<h3 class="date"><br>Session 5</h3>
<h2 class="author"><br>Dr. Richard M. Crowley<br>rcrowley@smu.edu.sg<br><a href="http://rmc.link/">http://rmc.link/</a></h2>

<!-- Load in primary data set -->
```{r, include=F, eval=F}
crsp <- read.csv("CCM_daily.csv")
crsp$date <- as.Date(crsp$date)
crsp <- crsp %>%
  arrange(gvkey, date) %>%
  mutate(ret = (cshoc - lag(cshoc))/lag(cshoc),
         mve = cshoc * prccd)

crsp <- crsp %>%
  group_by(gvkey) %>%
  mutate(ret.sd = rollapply(data=ret, width=253, FUN=sd, align="right", fill=NA, na.rm=T)) %>%
  ungroup()
crsp <- full_join(crsp, df_dates)
crsp <- crsp %>%
  arrange(gvkey, date)
crsp <- crsp %>% group_by(gvkey) %>% mutate(ret.sd1 = ifelse(is.na(ret.sd) & from_df == 1, lag(ret.sd,1), ret.sd)) %>% ungroup()
crsp <- crsp %>% filter(from_df == 1) %>% select(gvkey, date, ret, mve, ret.sd)
write.csv(crsp, "../../Data/Session_6-5.csv", row.names = F)
```

```{r, include=FALSE}
library(tidyverse)
library(plotly)
library(lubridate)
df <- read.csv("../../Data/Session_5-1.csv", stringsAsFactors=FALSE)
df_ratings <- read.csv("../../Data/Session_5-2.csv", stringsAsFactors=FALSE)
df_mve <- read.csv("../../Data/Session_5-3.csv", stringsAsFactors=FALSE)
df_rf <- read.csv("../../Data/Session_5-4.csv", stringsAsFactors=FALSE)
df_stock <- read.csv("../../Data/Session_5-5.csv", stringsAsFactors=FALSE)
```

# Front matter {data-background="../Backgrounds/section.jpg" class="header"}

## Learning objectives {data-background="../Backgrounds/default.jpg" class="default"}

<div class="clearfix">
<div class="left">
<img src="../Course_Diagrams/Map-5.png" width="500px">
</div>
<div class="right objectives">
- <b>Theory:</b>
    - Academic research
- <b>Application:</b>
    - Predicting bankruptcy over the next year for US manufacturing firms
        - Extend to credit downgrades
- <b>Methodology:</b>
    - Logistic regression
    - Models from academic research
</div>
<div class="clearfix">

## Datacamp {data-background="../Backgrounds/default.jpg" class="default"}

> - Explore on your own
> - No specific required class this week

# Academic research {data-background="../Backgrounds/section.jpg" class="header"}

## History of academic research in accounting {data-background="../Backgrounds/default.jpg" class="default"}

- Academic research in accounting, as it is today, began in the 1960s
    - What we call *Positive Accounting Theory*
        - Positive theory: understanding how the world works
- Prior to the 1960s, the focus was on Prescriptive theory
    - How the world should work
- Accounting research builds on work from many fields:
    - Economics
    - Finance
    - Psychology
    - Econometrics
    - Computer science (more recently)

## Types of academic research {data-background="../Backgrounds/default.jpg" class="default"}

- Theory
    - Pure economics proofs and simulation
- Experimental
    - Proper experimentation done on individuals
    - Can be psychology experiments or economic experiments
- Empirical/Archival
    - Data driven research
    - Based on the usage of historical data (i.e., archives)
    - Most likely to be easily co-optable by businesses and regulators

## Empirical/archival research in accounting {data-background="../Backgrounds/default.jpg" class="default present-only"}

<a href="https://rmc.link/Slides/SMU_IS_2021/Slides/Slides_post.html">Example</a>

## Who leverages accounting research {data-background="../Backgrounds/default.jpg" class="default"}

- Hedge funds
- Mutual funds
- Auditors
- Law firms
- Government entities like SG MAS and US SEC
- Exchanges like SGX

## Where can you find academic research {data-background="../Backgrounds/default.jpg" class="default"}

- The <a href="https://library.smu.edu.sg/">SMU library</a> has access to almost all high quality accounting research
- <a href="http://scholar.google.com/">Google Scholar</a> is a great site to discover research past and present
- <a href="https://www.ssrn.com/en/">SSRN</a> is the site to find cutting edge accounting and business research
    - <a href="https://papers.ssrn.com/sol3/topten/topTenResults.cfm?groupingId=204&netorjrnl=ntwk">List of top accounting papers on SSRN</a> (by downloads)

# Academic models: Altman Z-Score {data-background="../Backgrounds/section.jpg" class="header"}

## First: Why care about bankruptcy? {data-background="../Backgrounds/group.jpg" class="default"}

- Read this article: <a href="https://www.supplychaindive.com/news/carillion-bankruptcy-supply-chain-problem-common-sources/516567/">rmc.link/420class5-1</a>
    - "Carillion's liquidation reveals the dangers of shared sourcing"

> Based on this article, why do we care about bankruptcy risk for other firms?

## Where does the model come from? {data-background="../Backgrounds/default.jpg" class="default"}

<div class="clearfix">
<div class="left">
- Altman 1968, Journal of Finance
- A seminal paper in Finance cited over 15,000 times by other academic papers
</div>
<div class="right">
<a href="https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1540-6261.1968.tb00843.x"><img src="../Figures/Altman1968.png" width="350px"></a>
</div>
</div>


## What is the model about? {data-background="../Backgrounds/default.jpg" class="default"}

- The model was developed to identify firms that are likely to go bankrupt out of a pool of firms
- Focuses on using ratio analysis to determine such firms

## Model specification {data-background="../Backgrounds/default.jpg" class="default"}

$$
Z  = 1.2 x_1 + 1.4 x_2 + 3.3 x_3 + 0.6 x_4 + 0.999 x_5
$$

- $x_1$: Working capital to assets ratio
- $x_2$: Retained earnings to assets ratio
- $x_3$: EBIT to assets ratio
- $x_4$: Market value of equity to book value of liabilities
- $x_5$: Sales to total assets

> This looks like a linear regression without a constant

## How did the measure come to be? {data-background="../Backgrounds/case.jpg" class="default"}

- It actually isn't a linear regression
    - It is a clustering method called MDA (multiple discriminant analysis)
        - There are newer methods these days, such as SVM
- Used data from 1946 through 1965
    - 33 US manufacturing firms that went bankrupt, 33 that survived

More about this, from Altman himself in 2000: <a href="http://pages.stern.nyu.edu/~ealtman/Zscores.pdf">rmc.link/420class5-2</a>

- Read the section "Variable Selection" starting on page 8
    - Skim through $x_1$, $x_2$, $x_3$, $x_4$, and $x_5$

> How would these assumptions stand today?

## Rating Altman's assumptions {data-background="../Backgrounds/rain.jpg" class="default present-only"}

<iframe allowfullscreen frameborder="0" height="720px" mozallowfullscreen src="https://app.wooclap.com/events/COOSTC/questions/5f571fc47b4d5b0e165526d0" width="100%"></iframe>

## Rating Altman's assumptions {data-background="../Backgrounds/case.jpg" class="default present-only never"}

<iframe allowfullscreen frameborder="0" height="600px" mozallowfullscreen src="https://app.wooclap.com/events/HZRQJS/questions/5f604092c1fe45642c8c2523" width="100%"></iframe>

## Who uses Altman Z? {data-background="../Backgrounds/default.jpg" class="default"}

- Despite the model's simplicity and age, it is still in use
    - The simplicity of it plays a large part
- Frequently used by financial analysts

> <a href="https://news.google.com/search?q=%22altman+z%22">Recent news mentioning it</a>

# Application {data-background="../Backgrounds/section.jpg" class="header"}

## Main question {data-background="../Backgrounds/default.jpg" class="default"}

> Can we use bankruptcy models to predict supplier bankruptcies?

But first:

> Does the Altman Z-score [still] pick up bankruptcy?

## Question structure {data-background="../Backgrounds/question.jpg" class="default"}

> Is this a forecasting or forensics question?

<div class="fragment exclude">

- It has a time dimension like a forecasting question
- It has a feeling of a forensics question

</div>

## The data {data-background="../Backgrounds/default.jpg" class="default"}

- Compustat provides data on bankruptcies, including the date a company went bankrupt
    - Bankruptcy information is included in the "footnote" items in Compustat
        - If `dlsrn == 2`, then the firm went bankrupt
        - Bankruptcy date is `dldte`
- Most components of the Altman Z-Score model are in Compustat
    - But we'll pull market value from CRSP, since it is more complete
- All components of our later models are from Compustat as well
- Company credit rating data also from Compustat (Rankings)

## Bankruptcy in the US {data-background="../Backgrounds/default.jpg" class="default"}

- Chapter 7
    - The company ceases operating and liquidates
- Chapter 11
    - For firms intending to reorganize the company to "try to become profitable again" (<a href="https://www.sec.gov/reportspubs/investor-publications/investorpubsbankrupthtm.html">US SEC</a>)

<img src="../Figures/Bankrupt_Businessman_Cartoon.svg" height="200px">

## Common outcomes of bankruptcy {data-background="../Backgrounds/default.jpg" class="default"}

1. Cease operations entirely (liquidated)
    - In which case the assets are often sold off
2. Acquired by another company
3. Merge with another company
4. Successfully restructure and continue operating as the same firm
5. Restructure and operate as a new firm

<img src="../Figures/Bankrupt_Businessman_Cartoon.svg" height="200px">

## Calculating bankruptcy {data-background="../Backgrounds/coding.jpg" class="default"}

```{r}
# initial cleaning
# 100338 is an outlier in the bonds distribution
df <- df %>% filter(at >= 1, revt >= 1, gvkey != 100338)

## Merge in stock value
df$date <- as.Date(df$datadate)
df_mve <- df_mve %>%
  mutate(date = as.Date(datadate),
         mve = csho * prcc_f) %>%
  rename(gvkey=GVKEY)

df <- left_join(df, df_mve[,c("gvkey","date","mve")])

df <- df %>%
  group_by(gvkey) %>%
  arrange(datadate) %>%
  mutate(bankrupt = ifelse(row_number() == n() & dlrsn == 2 &
                           !is.na(dlrsn), 1, 0),
         bankrupt_lead = lead(bankrupt)) %>%
  ungroup() %>%
  filter(!is.na(bankrupt_lead)) %>%
  mutate(bankrupt_lead = factor(bankrupt_lead, levels=c(0,1)))
```

- `row_number()` gives the current row within the group, with the first row as 1, next as 2, etc.
- `n()` gives the number of rows in the group

## Calculating the Altman Z-Score {data-background="../Backgrounds/coding.jpg" class="default"}

```{r}

# Calculate the measures needed
df <- df %>%
  mutate(wcap_at = wcap / at,  # x1
         re_at = re / at,  # x2
         ebit_at = ebit / at,  # x3
         mve_lt = mve / lt,  # x4
         revt_at = revt / at)  # x5
# cleanup
df <- df %>%
  mutate_if(is.numeric, list(~replace(., !is.finite(.), NA)))

# Calculate the score
df <- df %>%
  mutate(Z = 1.2 * wcap_at + 1.4 * re_at + 3.3 * ebit_at + 0.6 * mve_lt + 
           0.999 * revt_at)

# Calculate date info for merging
df$date <- as.Date(df$datadate)
df$year <- year(df$date)
df$month <- month(df$date)
```

- Calculate $x_1$ through $x_5$
- Apply the model directly

## Build in credit ratings {data-background="../Backgrounds/coding.jpg" class="default"}

> We'll check our Z-score against credit rating as a simple validation

```{r}
# df_ratings has ratings data in it

# Ratings, in order from worst to best
ratings <- c("D", "C", "CC", "CCC-", "CCC","CCC+", "B-", "B", "B+", "BB-",
             "BB", "BB+", "BBB-", "BBB", "BBB+", "A-", "A", "A+", "AA-", "AA",
             "AA+", "AAA-", "AAA", "AAA+")
# Convert string ratings (splticrm) to numeric ratings
df_ratings$rating <- factor(df_ratings$splticrm, levels=ratings, ordered=T)

df_ratings$date <- as.Date(df_ratings$datadate)
df_ratings$year <- year(df_ratings$date)
df_ratings$month <- month(df_ratings$date)

# Merge together data
df <- left_join(df, df_ratings[,c("gvkey", "year", "month", "rating")])
```

## Z vs credit ratings, 1973-2017 {data-background="../Backgrounds/case.jpg" class="default"}

<div class="clearfix">
<div class="right">
```{r, echo=F, fig.height=5, fig.width=6}
plot <- df %>%
  filter(!is.na(Z), !is.na(rating)) %>%
  group_by(rating) %>%
  mutate(mean_Z=mean(Z,na.rm=T)) %>%
  slice(1) %>%
  ungroup() %>%
  select(rating, mean_Z) %>%
  ggplot(aes(y=mean_Z, x=rating)) + 
  geom_col() + 
  ylab('Mean Altman Z') + xlab('Credit rating') + 
  theme(axis.text.x = element_text(angle = 90))
ggplotly(plot)
```
</div>
<div class="right">
```{r}
df %>%
  filter(!is.na(Z),
         !is.na(bankrupt)) %>%
  group_by(bankrupt_lead) %>%
  mutate(mean_Z=mean(Z,na.rm=T)) %>%
  slice(1) %>%
  ungroup() %>%
  select(bankrupt_lead, mean_Z) %>%
  html_df()
```
</div>

## Z vs credit ratings, 2000-2017 {data-background="../Backgrounds/case.jpg" class="default"}

<div class="clearfix">
<div class="right">
```{r, echo=F, fig.height=5, fig.width=6}
plot <- df %>%
  filter(!is.na(Z), !is.na(rating), year >= 2000) %>%
  group_by(rating) %>%
  mutate(mean_Z=mean(Z,na.rm=T)) %>%
  slice(1) %>%
  ungroup() %>%
  select(rating, mean_Z) %>%
  ggplot(aes(y=mean_Z, x=rating)) + 
  geom_col() + 
  ylab('Mean Altman Z') + xlab('Credit rating') + 
  theme(axis.text.x = element_text(angle = 90))
ggplotly(plot)
```
</div>
<div class="right">
```{r}
df %>%
  filter(!is.na(Z),
         !is.na(bankrupt_lead),
         year >= 2000) %>%
  group_by(bankrupt_lead) %>%
  mutate(mean_Z=mean(Z,na.rm=T)) %>%
  slice(1) %>%
  ungroup() %>%
  select(bankrupt_lead, mean_Z) %>%
  html_df()
```
</div>

## Test it with a regression {data-background="../Backgrounds/case.jpg" class="default"}

```{r, warning=F}
fit_Z <- glm(bankrupt_lead ~ Z, data=df, family=binomial)
summary(fit_Z)
```

## How good is the model though??? {data-background="../Backgrounds/default.jpg" class="default"}

Examples:

> Correct 92.6% of the time using Z < 1 as a cutoff

- Only correctly captures 29 of 85 bankruptcies

> Correct 99.8% of the time if we say firms never go bankrupt...

- Correctly captures **0** of 85 bankruptcies

```{r, echo=F}
lz <- df %>% filter(!is.na(bankrupt_lead), !is.na(Z)) %>% filter(Z < 1) %>% pull(bankrupt_lead) %>% table()
hz <- df %>% filter(!is.na(bankrupt_lead), !is.na(Z)) %>% filter(Z >= 1) %>% pull(bankrupt_lead) %>% table()
x <- matrix(c(lz, hz), nrow=2)
rownames(x) <- c('No bankruptcy', 'Bankruptcy')
colnames(x) <- c('Z < 1', 'Z >= 1')
x
```

# Errors in binary testing {data-background="../Backgrounds/section.jpg" class="header"}

## Types of errors {data-background="../Backgrounds/reference.jpg" class="default"}

<img src="../Figures/error_types.png" height="450px">

> This type of chart (filled in) is called a *Confusion matrix*

## Type I error (False positive) {data-background="../Backgrounds/reference.jpg" class="default"}

> We say that the company will go bankrupt, but they don't

- A Type I error occurs any time we say something is *true*, yet it is false
- Quantifying type I errors in the data
    - False positive rate (FPR)
        - The percent of failures misclassified as successes
    - Specificity: $1 - FPR$
        - A.k.a. true negative rate (TNR)
        - The percent of failures properly classified

<img src="../Figures/FPR.png">

## Type 2 error (False negative) {data-background="../Backgrounds/reference.jpg" class="default"}

> We say that the company *will not* go bankrupt, yet they do

- A Type II error occurs any time we say something is *false*, yet it is true
- Quantifying type I errors in the data
    - False negative rate (FNR): $1-Sensitivity$
        - The percent of successes misclassified as failures
    - Sensitivity:
        - A.k.a. true positive rate (TPR)
        - The percent of successes properly classified

<img src="../Figures/TPR.png">

## Useful equations {data-background="../Backgrounds/reference.jpg" class="default"}


<img src="../Figures/FPR.png" width="600px">
<img src="../Figures/TPR.png" width="600px">
<img src="../Figures/TNR.png" width="600px">
<img src="../Figures/Accuracy.png" width="600px">

## A note on the equations {data-background="../Backgrounds/reference.jpg" class="default"}

- Accuracy is very useful if you are predicting something that occurs reasonably frequently
    - Not too often, but not too rarely
- Sensitivity is very useful for rare events
- Specificity is very useful for frequent events
    - Or for events where misclassifying the null is very troublesome
        - Criminal trials
        - Medical diagnoses

> Calculating any of these require the following

1. Predict percentages using `predict(. , type="response")`
2. Convert predictions to binary outcomes by specifying a desired cutoff
    - Anything below the cutoff is `0`, anything above it is `1`

## A more comprehensive approach {data-background="../Backgrounds/case_digital.jpg" class="default"}

- Using `package:yardstick` we can plot out specificity and sensitivity across all possible cutoffs!

```{r, message=F, error=F, warning=F, fig.height=5, fig.width=5}
library(yardstick)
df_Z <- df %>% filter(!is.na(Z), !is.na(bankrupt_lead))
df_Z$pred <- predict(fit_Z, df_Z, type="response")
df_Z %>% roc_curve(truth=bankrupt_lead, estimate=pred, event_level='second') %>%
  autoplot()
```



## ROC curves {data-background="../Backgrounds/reference.jpg" class="default"}

- The previous graph is called a ROC curve, or *r*eceiver *o*perator *c*haracteristic curve
- The higher up and left the curve is, the better the logistic regression fits.

<div class="clearfix">
<div class="left">
- Neat properties:
    - The area under a perfect model is always 1
    - The area under random chance is always 0.5
        - This is the straight dashed line on the graph
</div>
<div class="right">
```{r, echo=F, fig.width=5, fig.height=5}
df_Z %>% roc_curve(truth=bankrupt_lead,
                  estimate=pred,
                  event_level='second') %>%
  autoplot()
```
</div>
</div>

## ROC AUC {data-background="../Backgrounds/reference.jpg" class="default"}

- The neat properties of the curve give rise to a useful statistic: ROC AUC
    - AUC = Area under the curve
- Ranges from 0 (perfectly incorrect) to 1 (perfectly correct)
- Above 0.6 is generally the minimum acceptable bound
    - 0.7 is preferred
    - 0.8 is very good
- `package:yardstick` can calculate this too

```{r}
auc_Z <- df_Z %>% roc_auc(truth=bankrupt_lead, estimate=pred, event_level='second')
auc_Z
```

## ROC AUC simplest interpretation {data-background="../Backgrounds/reference.jpg" class="default"}

<img src="../Figures/AUC_Interpretation.png" width="600px">

> AUC is the probability that our model assigns a higher estimated probability to a randomly selected 1 than to a randomly selected 0.

## R Practice ROC AUC {data-background="../Backgrounds/coding.jpg" class="default"}

- Practice using these new functions with last week's Walmart data
    1. Model decreases in revenue using prior quarter YoY revenue growth
    2. Explore the model using `predict()`
    3. Calculate ROC AUC
    4. Plot a ROC curve
- Do all exercises in today's practice file
    - <a href="http://rmc.link/Slides/acct420v4/Session_5/Session_5_R.html">R Practice</a>
    - Shortlink: <a href="http://rmc.link/Slides/acct420v4/Session_5/Session_5_R.html">rmc.link/420r5</a>

# Academic models: Distance to default (DD) {data-background="../Backgrounds/section.jpg" class="header"}

## Where does the model come from? {data-background="../Backgrounds/default.jpg" class="default"}

<div class="clearfix">
<div class="left">
- Merton 1974, Journal of Finance
- Another seminal paper in finance, cited by over 12,000 other academic papers
- <a href="https://www.nobelprize.org/prizes/economics/1997/merton/facts/">About Merton</a>
</div>
<div class="right">
<a href="https://onlinelibrary.wiley.com/doi/full/10.1111/j.1540-6261.1974.tb03058.x"><img src="../Figures/merton1974.png" width="400px"></a>
</div>
</div>

## What is the model about? {data-background="../Backgrounds/default.jpg" class="default"}

- The model itself comes from thinking of debt in an options pricing framework
- Uses the Black-Scholes model to price out a company
- Consider a company to be bankrupt when the company is not worth more than the the debt itself, in expectation

## Model specification {data-background="../Backgrounds/default.jpg" class="default"}

$$
DD = \frac{\log(V_A / D) + (r-\frac{1}{2}\sigma_A^2)(T-t)}{\sigma_A \sqrt(T-t)}
$$

<div class="clearfix">
<div class="left">
- $V_A$: Value of assets
    - Market based
- $D$: Value of liabilities
    - From balance sheet
- $r$: The risk free rate
- $\sigma_A$: Volatility of assets
    - Use daily stock return volatility, annualized
        - Annualized means multiply by $\sqrt{253}$
- $T-t$: Time horizon
</div>
<div class="right">
```{r, echo=F}
score = 1
m = 0
std = 1

funcShaded <- function(x, lower_bound) {
    y = dnorm(x, mean = m, sd = std)
    y[x < lower_bound] <- NA
    return(y)
}

ggplot(data.frame(x = c(-3, 3)), aes(x = x)) + 
  stat_function(fun = dnorm, args = list(mean = m, sd = std)) + 
  stat_function(fun = funcShaded, args = list(lower_bound = score), 
                geom = "area", fill = 'black', alpha = .2) +
  scale_x_continuous(name = "Score", breaks = seq(-3, 3, std)) + 
  geom_text(data = data.frame(x=c(1.5), y=c(0.05)), aes(x=x, y=y, label="Prob(default)", size=30)) + 
  geom_line(data = data.frame(x=c(1,1), y=c(0,0.4)), aes(x=x,y=y)) + 
  geom_text(data = data.frame(x=c(1.3), y=c(0.4)), aes(x=x, y=y, label="DD", size=30)) +
  theme(legend.position="none")
```
</div>
</div>


## Who uses it? {data-background="../Backgrounds/default.jpg" class="default"}

- Moody's KMV is derived from the Merton model
    - Common platform for analyzing risk in financial services
    - <a href="https://www.moodysanalytics.com/solutions-overview/credit-risk/credit-risk-modeling">More information</a>

<img src="../Figures/moodys.png" height="200px">

# Applying DD {data-background="../Backgrounds/section.jpg" class="header"}

## Calculating DD in R {data-background="../Backgrounds/coding.jpg" class="default"}

- First we need one more measure: the standard deviation of assets
    - This varies by time, and construction of it is subjective
    - We will use standard deviation over the last 5 years

```{r}
# df_stock is an already prepped csv from CRSP data
df_stock$date <- as.Date(df_stock$date)
df <- left_join(df, df_stock[,c("gvkey", "date", "ret", "ret.sd")])
```

## Calculating DD in R {data-background="../Backgrounds/coding.jpg" class="default"}

```{r}

df_rf$date <- as.Date(df_rf$dateff)
df_rf$year <- year(df_rf$date)
df_rf$month <- month(df_rf$date)

df <- left_join(df, df_rf[,c("year", "month", "rf")])

df <- df %>%
  mutate(DD = (log(mve / lt) + (rf - (ret.sd*sqrt(253))^2 / 2)) /
              (ret.sd*sqrt(253)))
# Clean the measure
df <- df %>%
  mutate_if(is.numeric, list(~replace(., !is.finite(.), NA)))
```

- Just apply the formula using mutate
- $\sqrt{253}$ is included because `ret.sd` is daily return standard deviation
    - There are ~253 trading days per year in the US

## DD vs credit ratings, 1973-2017 {data-background="../Backgrounds/case.jpg" class="default"}

<div class="clearfix">
<div class="right">
```{r, echo=F, fig.height=5, fig.width=4}
plot <- df %>%
  filter(!is.na(DD),
         !is.na(rating)) %>%
  group_by(rating) %>%
  mutate(mean_DD=mean(DD,na.rm=T),
         prob_default = pnorm(-1 * mean_DD)) %>%
  slice(1) %>%
  ungroup() %>%
  select(rating, prob_default) %>%
  ggplot(aes(y=prob_default, x=rating)) + 
  geom_col() + 
  ylab('Probability of default') + xlab('Credit rating') + 
  theme(axis.text.x = element_text(angle = 90))
ggplotly(plot)
```
</div>
<div class="left">
```{r}
df %>%
  filter(!is.na(DD),
         !is.na(bankrupt_lead)) %>%
  group_by(bankrupt_lead) %>%
  mutate(mean_DD=mean(DD, na.rm=T),
         prob_default =
           pnorm(-1 * mean_DD)) %>%
  slice(1) %>%
  ungroup() %>%
  select(bankrupt_lead, mean_DD,
         prob_default) %>%
  html_df()
```
</div>

## DD vs credit ratings, 2000-2017 {data-background="../Backgrounds/case.jpg" class="default"}

<div class="clearfix">
<div class="right">
```{r, echo=F, fig.height=5, fig.width=4}
plot <- df %>%
  filter(!is.na(DD),
         !is.na(rating),
         year >= 2000) %>%
  group_by(rating) %>%
  mutate(mean_DD=mean(DD,na.rm=T),
         prob_default = pnorm(-1 * mean_DD)) %>%
  slice(1) %>%
  ungroup() %>%
  select(rating, prob_default) %>%
  ggplot(aes(y=prob_default, x=rating)) + 
  geom_col() + 
  ylab('Probability of default') + xlab('Credit rating') + 
  theme(axis.text.x = element_text(angle = 90))
ggplotly(plot)
```
</div>
<div class="left">
```{r}
df %>%
  filter(!is.na(DD),
         !is.na(bankrupt_lead),
         year >= 2000) %>%
  group_by(bankrupt_lead) %>%
  mutate(mean_DD=mean(DD, na.rm=T),
         prob_default =
           pnorm(-1 * mean_DD)) %>%
  slice(1) %>%
  ungroup() %>%
  select(bankrupt_lead, mean_DD,
         prob_default) %>%
  html_df()
```
</div>

## Test it with a regression {data-background="../Backgrounds/case.jpg" class="default"}

```{r, warning=FALSE}
fit_DD <- glm(bankrupt_lead ~ DD, data=df, family=binomial)
summary(fit_DD)
```

## ROC Curves {data-background="../Backgrounds/case.jpg" class="default"}

```{r, fig.height=3.5, fig.width=3.5}
df_DD <- df %>% filter(!is.na(Z), !is.na(bankrupt_lead))
df_DD$pred <- predict(fit_DD, df_DD, type="response")
df_DD %>% roc_curve(truth=bankrupt_lead, estimate=pred, event_level='second') %>%
  autoplot()
df_DD %>% roc_auc(truth=bankrupt_lead, estimate=pred, event_level='second')
```

## AUC comparison {data-background="../Backgrounds/case.jpg" class="default"}

```{r}
#AUC
auc_DD <- df_DD %>% roc_auc(truth=bankrupt_lead, estimate=pred, event_level='second')
AUCs <- c(auc_Z$.estimate, auc_DD$.estimate)
names(AUCs) <- c("Z", "DD")
AUCs
```

> Distance to Default performs quite a bit better here!

# A more practical application {data-background="../Backgrounds/section.jpg" class="header"}

## A more practical application {data-background="../Backgrounds/default.jpg" class="default"}

- Companies don't only have problems when there is a bankruptcy
    - Credit downgrades can be just as bad

> Why?

<div class="exclude">

- Credit downgrades cause an increase in interest rates for debt, leading to potential liquidity issues.

</div>

## Predicting downgrades {data-background="../Backgrounds/coding.jpg" class="default"}

```{r, warning=FALSE}
# calculate downgrade
df <- df %>%
  group_by(gvkey) %>%
  arrange(date) %>%
  mutate(downgrade = factor(ifelse(lead(rating) < rating,1,0), levels=c(0,1)),
         diff_Z = Z - lag(Z),
         diff_DD = DD - lag(DD)) %>%
  ungroup()


# training sample
train <- df %>% filter(year < 2014, !is.na(diff_Z), !is.na(diff_DD), !is.na(downgrade),
                       year > 1985)
test <- df %>% filter(year >= 2014, !is.na(diff_Z), !is.na(diff_DD), !is.na(downgrade))

# glms
fit_Z2 <- glm(downgrade ~ diff_Z, data=train, family=binomial)
fit_DD2 <- glm(downgrade ~ diff_DD, data=train, family=binomial)
```

## Predicting downgrades with Altman Z {data-background="../Backgrounds/case.jpg" class="default"}

```{r}
summary(fit_Z2)
```


## Predicting downgrades with DD {data-background="../Backgrounds/case.jpg" class="default"}

```{r}
summary(fit_DD2)
```

## ROC Performance on this task {data-background="../Backgrounds/case.jpg" class="default"}

```{r, echo=F, fig.height=6, fig.width=8}
df_Z2 <- train %>% filter(!is.na(diff_Z), !is.na(downgrade))
df_Z2$pred <- predict(fit_Z2, df_Z2, type="response")
curve1 <- df_Z2 %>% roc_curve(truth=downgrade, estimate=pred, event_level='second')
auc_Z2 <- df_Z2 %>% roc_auc(truth=downgrade, estimate=pred, event_level='second')

df_DD2 <- train %>% filter(!is.na(diff_DD), !is.na(downgrade))
df_DD2$pred <- predict(fit_DD2, df_DD2, type="response")
curve2 <- df_DD2 %>% roc_curve(truth=downgrade, estimate=pred, event_level='second')
auc_DD2 <- df_DD2 %>% roc_auc(truth=downgrade, estimate=pred, event_level='second')

curve1 <- curve1 %>% group_by(sensitivity) %>% slice(c(1, n())) %>% ungroup()
curve2 <- curve2 %>% group_by(sensitivity) %>% slice(c(1, n())) %>% ungroup()

ggplot() +
  geom_line(data=curve1, aes(y=sensitivity, x=1-specificity, color="Altman Z")) + 
  geom_line(data=curve2, aes(y=sensitivity, x=1-specificity, color="DD")) + 
  geom_abline(slope=1)
```



```{r, echo=F}
AUCs <- c(auc_Z2$.estimate, auc_DD2$.estimate)
names(AUCs) <- c("Z", "DD")
AUCs
```

## Out of sample ROC performance {data-background="../Backgrounds/case.jpg" class="default"}

```{r, echo=F, fig.height=5}
df_Z2 <- test %>% filter(!is.na(diff_Z), !is.na(downgrade))
df_Z2$pred <- predict(fit_Z2, df_Z2, type="response")
curve1 <- df_Z2 %>% roc_curve(truth=downgrade, estimate=pred, event_level='second')
auc_Z2 <- df_Z2 %>% roc_auc(truth=downgrade, estimate=pred, event_level='second')

df_DD2 <- test %>% filter(!is.na(diff_DD), !is.na(downgrade))
df_DD2$pred <- predict(fit_DD2, df_DD2, type="response")
curve2 <- df_DD2 %>% roc_curve(truth=downgrade, estimate=pred, event_level='second')
auc_DD2 <- df_DD2 %>% roc_auc(truth=downgrade, estimate=pred, event_level='second')

curve1 <- curve1 %>% group_by(sensitivity) %>% slice(c(1, n())) %>% ungroup()
curve2 <- curve2 %>% group_by(sensitivity) %>% slice(c(1, n())) %>% ungroup()

ggplot() +
  geom_line(data=curve1, aes(y=sensitivity, x=1-specificity, color="Altman Z")) + 
  geom_line(data=curve2, aes(y=sensitivity, x=1-specificity, color="DD")) + 
  geom_abline(slope=1)
```

```{r, echo=F}
AUCs <- c(auc_Z2$.estimate, auc_DD2$.estimate)
names(AUCs) <- c("Z", "DD")
AUCs
```

## Predicting bankruptcy {data-background="../Backgrounds/group.jpg" class="default"}

> What other data could we use to predict corporate bankruptcy as it relates to a company's supply chain?

- What is the reason that this event or data would be useful for prediction?
    - I.e., how does it fit into your mental model?
- A useful starting point from McKinsey
    - <a href="https://www.mckinsey.com/business-functions/operations/our-insights/big-data-and-the-supply-chain-the-big-supply-chain-analytics-landscape-part-1">rmc.link/420class5-3</a>
        - Section "B. Sourcing"

## Brainstorming... {data-background="../Backgrounds/rain.jpg" class="default present-only"}

<iframe allowfullscreen frameborder="0" height="600px" mozallowfullscreen src="https://app.wooclap.com/events/HZRQJS/questions/5f605a7acd510a640f109fc1" width="100%"></iframe>

# Combined model {data-background="../Backgrounds/section.jpg" class="header"}

## Building a combined model {data-background="../Backgrounds/default.jpg" class="default"}

```{r, warning=F}
fit_comb <- glm(downgrade ~ diff_Z + diff_DD, data=train, family=binomial)
summary(fit_comb)
```

## Marginal effects {data-background="../Backgrounds/default.jpg" class="default"}

```{r}
library(margins)
fit_comb %>%
  margins() %>%
  summary()
```

## ROC Performance on this task {data-background="../Backgrounds/case.jpg" class="default"}

```{r, echo=F, fig.height=5}
df_comb <- test %>% filter(!is.na(diff_DD), !is.na(diff_Z), !is.na(downgrade))
df_comb$pred <- predict(fit_comb, df_comb, type="response")
curve3 <- df_comb %>% roc_curve(truth=downgrade, estimate=pred, event_level='second')
auc_comb <- df_comb %>% roc_auc(truth=downgrade, estimate=pred, event_level='second')

curve3 <- curve3 %>% group_by(sensitivity) %>% slice(c(1, n())) %>% ungroup()

ggplot() +
  geom_line(data=curve1, aes(y=sensitivity, x=1-specificity, color="Altman Z")) + 
  geom_line(data=curve2, aes(y=sensitivity, x=1-specificity, color="DD")) +
  geom_line(data=curve3, aes(y=sensitivity, x=1-specificity, color="Combined")) + 
  geom_abline(slope=1)

AUCs <- c(auc_Z2$.estimate, auc_DD2$.estimate, auc_comb$.estimate)
names(AUCs) <- c("Z", "DD", "Combined")
AUCs
```

# End matter {data-background="../Backgrounds/ending.jpg" class="header"}

## For next week {data-background="../Backgrounds/ending.jpg" class="default"}

- For next week:
    - Second individual assignment
        - Finish by the end of the day next class session
        - Submit on eLearn
    - Datacamp
        - Practice a bit more to keep up to date
            - Using R more will make it more natural
    - Keep thinking about who you want to work with on the project
- Survey on the class session at <a href="https://docs.google.com/forms/d/e/1FAIpQLSeCqZUD2NeYkavacSj6T4RVMlZMWiujDG_M7lmvL8esr2GPXg/viewform?usp=sf_link">rmc.link/420survey5</a>

## Packages used for these slides {data-background="../Backgrounds/code.jpg" class="default"}

- `package:kableExtra`
- `package:knitr`
- `package:lubridate`
- `package:magrittr`
- `package:margins`
- `package:plotly`
- `package:revealjs`
- `package:tidyverse`
- `package:yardstick`

## References {data-background="../Backgrounds/books-old.jpg" class="default stack never"}

## Custom code {data-background="../Backgrounds/code.jpg" class="default"}

<small style="width:100%">
```{r, eval=F}
# Calculating a 253 day rolling standard deviation in R
crsp <- crsp %>%
  group_by(gvkey) %>%
  mutate(ret.sd = rollapply(data=ret, width=253, FUN=sd, align="right", fill=NA, na.rm=T)) %>%
  ungroup()
```
</small>


## NOTES {class="never"}

