---
title: "ACCT 420: Logistic Regression for Bankruptcy"
author: 
  - name: Dr. Richard M. Crowley
    email: rcrowley@smu.edu.sg
    affiliation: https://rmc.link/
title-slide-attributes:
  data-background-video: "../Backgrounds/title.mp4"
  data-background-video-loop: "true"
  data-background-size: stretch
format:
  revealjs:
    theme: default
    css: quarto_white.css
    navigation-mode: vertical
    controls: true
    transition: slide
    background-transition: fade
    transition-speed: default
    width: 1536
    height: 864
    logo: ../Logos/logo2.png
    code-link: true
    code-line-numbers: false
    chalkboard:
      theme: whiteboard
    auto-stretch: false
revealjs-plugins:
  - attribution
execute:
  echo: true
  warning: false
  message: false
editor: source
---

# Front Matter {data-background-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

<!--# Note that the downlit and xml2 pacakges are required to compile -->
<!--# Installed extensions include `attribution`, `qrcode`, and `material-icons` -->

```{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_5-5.csv", row.names = F)
```

```{r helpers, echo=FALSE, warning=FALSE, message=F}
library(knitr)
library(kableExtra)

#' Pretty formatted HTML tables
#' 
#' @param df The dataframe to format as a table
#' @param names Names to use for columns if not using the df column names. Vector
#' @param highlight_cols Columns to highlight. Default is c(1) for highlighting the first row. Use c() to remove.
#' @param highlight_rows Rows to highlight. Default is c() to not highlight any rows.
#' @param color_rows Rows to highlight in color. Default is c() to not highlight any rows.
#' @param highlight_end_rows Rows to highlight. Uses a reverse indexing where 1 is the last row, 2 is second to last, etc.
#' @param full Use to span the full width of a slide. FALSE by default.
#' @param fixed_header Whether the first row should be treated as a header row. TRUE by default.
html_df <- function(df, names=NULL, highlight_cols=c(1), highlight_rows=c(), color_rows=c(), highlight_end_rows=c(), full=F, fixed_header=T) {
  if(!length(names)) {
    names=colnames(df)
  }
  kbl(df,"html", col.names = names, align = c("l",rep('c',length(cols)-1))) %>%
    kable_paper(c("striped","hover"), full_width=full, fixed_thead=F, html_font = "\"Source Sans Pro\", Helvetica, sans-serif") %>%
    {if(length(highlight_cols)) column_spec(., highlight_cols,bold=T, background = "#ffffff99") else .} %>%
    {if(length(highlight_rows)) row_spec(., highlight_rows,bold=T, background = "#ffffff99") else .} %>%
    {if(length(color_rows)) row_spec(., color_rows,bold=T, background = "#aaaaff99") else .} %>%
    {if(length(highlight_end_rows)) row_spec(., nrow(df) + 1 - highlight_end_rows,bold=T, background = "#ffffff99") else .}
}

library(tidyverse)
```

```{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)
```

<!--# List packages here that are not being correctly highlighted -->
<!--# The usual cause is calling the library from a echo=FALSE cell -->

```{r}
#| include: false
options(downlit.attached = c("knitr", "kableExtra", "tidyverse"))
```

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

:::: columns

::: {.column width="50%"}

![](../Course_Diagrams/Map-5.png){fig-alt="Roadmap for the course"}

:::

::: {.column width="50%"}

- **Theory:**
    - Academic research
- **Application:**
    - Predicting bankruptcy over the next year for US manufacturing firms
        - Extend to credit downgrades
- **Methodology:**
    - Logistic regression
    - Models from academic research

:::

::::


# Academic research {data-background-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

- 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"}

- 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"}

[Click for examples](https://rmc.link/Slides/SMU_PhD_Workshop_2022/Reveal.js/Slides/Slides_post.html)

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


- 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"}

- 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-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

- 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?

## 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/613ed3c77c07621214fba1c3" width="100%"></iframe>


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

:::: {.columns}

::: {.column width="50%"}

- Altman 1968, Journal of Finance
- A seminal paper in Finance cited over 15,000 times by other academic papers

:::

::: {.column width="50%"}

[![](../Figures/Altman1968.png){width="350px" fig-align="center"}](https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1540-6261.1968.tb00843.x)

:::

::::

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

- 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"}

$$
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"}

- 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="600px" mozallowfullscreen src="https://app.wooclap.com/events/HZRQJS/questions/5f604092c1fe45642c8c2523" width="100%"></iframe>

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

- 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-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

> 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"}

> 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"}

- 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"}

- 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>)

![](../Figures/Bankrupt_Businessman_Cartoon.svg){height="200px" fig-align="center"}

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

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

![](../Figures/Bankrupt_Businessman_Cartoon.svg){height="200px" fig-align="center"}

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

```{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 <- ymd(df$datadate)
df_mve <- df_mve %>%
  mutate(date = ymd(datadate),
         mve = csho * prcc_f) %>%
  rename(gvkey=GVKEY)

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

```{r}
# Compute the bankruptcy indicator
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; the first row is 1, next is 2, etc.
- `n()` gives the number of rows in the group

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

```{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(across(where(is.numeric), ~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 <- ymd(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"}

> 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")])
```

> Credit rating is from Compustat -- this data is no longer provided

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

:::: {.columns}

::: {.column width="50%"}

```{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)
```

:::

::: {.column width="50%"}

```{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()
```

:::

::::

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

:::: {.columns}

::: {.column width="50%"}

```{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)
```

:::

::: {.column width="50%"}

```{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()
```

:::

::::

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

```{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"}

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-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

![](../Figures/error_types.png){height="450px" fig-align="center"}

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

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

> 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

![](../Figures/FPR.png){width="1200px" fig-align="center"}

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

> 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

![](../Figures/TPR.png){width="1200px" fig-align="center"}

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

![](../Figures/FPR.png){width="800px" fig-align="center"}
![](../Figures/TPR.png){width="800px" fig-align="center"}
![](../Figures/TNR.png){width="800px" fig-align="center"}
![](../Figures/Accuracy.png){width="800px" fig-align="center"}

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

- 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"}

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

```{r, message=F, error=F, warning=F, fig.height=4.8, fig.width=4.8}
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(bankrupt_lead, pred, event_level='second') %>% autoplot()
```



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

- 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.

:::: {.columns}

::: {.column width="50%"}

- 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

:::

::: {.column width="50%"}

```{r, echo=F, fig.width=4.8, fig.height=4.8}
df_Z %>% roc_curve(bankrupt_lead,
                  pred,
                  event_level='second') %>%
  autoplot()
```

:::

::::


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

- 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
- `{yardstick}` can calculate this too

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

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

![](../Figures/AUC_Interpretation.png){width="600px" fig-align="center"}

> 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"}

- 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/420r5">R Practice</a>
    - Shortlink: <a href="http://rmc.link/420r5">rmc.link/420r5</a>

# Academic models: Distance to default (DD) {data-background-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

:::: {.columns}

::: {.column width="50%"}

- Merton 1974, Journal of Finance
- Another seminal paper in finance, cited by over 12,000 other academic papers
- [About Merton](https://www.nobelprize.org/prizes/economic-sciences/1997/merton/facts/)

:::

::: {.column width="50%"}

[![](../Figures/merton1974.png){width="400px" fig-align="center"}](https://onlinelibrary.wiley.com/doi/full/10.1111/j.1540-6261.1974.tb03058.x)

:::

::::


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

- 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"}

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

:::: {.columns}

::: {.column width="50%"}

- $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

:::

::: {.column width="50%"}

```{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")
```

:::

::::

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

- 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>

![](../Figures/moodys.png){height="200px" fig-align="center"}

# Applying DD {data-background-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

- 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"}

```{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(across(where(is.numeric), ~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" .small}

:::: {.columns}

::: {.column width="50%"}

```{r, echo=F, fig.height=5, fig.width=6}
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)
```

:::

::: {.column width="50%"}

```{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()
```

:::

::::

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

:::: {.columns}

::: {.column width="50%"}

```{r, echo=F, fig.height=5, fig.width=6}
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)
```

:::

::: {.column width="50%"}

```{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()
```

:::

::::

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

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

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

```{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(bankrupt_lead, pred, event_level='second') %>%
  autoplot()
df_DD %>% roc_auc(bankrupt_lead, pred, event_level='second')
```

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

```{r}
#AUC
auc_DD <- df_DD %>% roc_auc(bankrupt_lead, 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-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

- 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"}

```{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"}

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


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

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

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

```{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(downgrade, pred, event_level='second')
auc_Z2 <- df_Z2 %>% roc_auc(downgrade, 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(downgrade, pred, event_level='second')
auc_DD2 <- df_DD2 %>% roc_auc(downgrade, 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"}

```{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(downgrade, pred, event_level='second')
auc_Z2 <- df_Z2 %>% roc_auc(downgrade, 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(downgrade, pred, event_level='second')
auc_DD2 <- df_DD2 %>% roc_auc(downgrade, 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"}

> 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
    - [rmc.link/420class5-3](https://www.mckinsey.com/business-functions/operations/our-insights/big-data-and-the-supply-chain-the-big-supply-chain-analytics-landscape-part-1)
        - 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-video="../Backgrounds/section.mp4" data-background-video-loop="true"}

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

```{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"}

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

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

```{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(downgrade, pred, event_level='second')
auc_comb <- df_comb %>% roc_auc(downgrade, 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-video="../Backgrounds/section.mp4"}

## Wrap up {data-background="../Backgrounds/ending.jpg"}

- 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 this QR code:

{{< qrcode https://rmc.link/420survey5 qrsurvey width=200 height=200 >}}

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

- `{DT}`
- `{downlit}`
- `{kableExtra}`
- `{knitr}`
- `{margins}`
- `{plotly}`
- `{quarto}`
- `{revealjs}`
- `{tidyverse}`
- `{yardstick}`


## Custom code {data-background="../Backgrounds/code.jpg" class="small" visibility="uncounted"}

```{r}
#| eval: false
#| echo: true

# 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()
```

## NOTES {class="never"}

- Replace the WRDS data using "Capital IQ - Security Ratings"
    - May require a good bit of changes to data

