---
title: "skills overlap"
author: "Jenna"
date: "June 17, 2020"
output: rmarkdown::github_document
---

```{r setup, include=FALSE}
#knitr::opts_chunk$set(echo = TRUE)
knitr::opts_chunk$set(echo = FALSE, warning=FALSE, message=FALSE)

```

```{r}
# Load 
library(tidyverse)
library(blscrapeR)
library(blsAPI)
library(fredr)
library(kableExtra)

library(here)
source(here("03_graphics themes.R"))

```

```{r}
blsOImatrix <- read_csv("inputs/skills/bls OI matrix.csv")
skills <- read_csv("inputs/skills/Skills.csv")
food <- read_csv("inputs/skills/food.csv") %>% mutate(NAICS = "72") #%>% top_n(n = 10, wt = `Employed in Industry (2018)`)
arts <- read_csv("inputs/skills/arts.csv") %>% mutate(NAICS = "71") #%>% top_n(n = 10, wt = `Employed in Industry (2018)`)
other <- read_csv("inputs/skills/other.csv") %>% mutate(NAICS = "81") #%>% top_n(n = 10, wt = `Employed in Industry (2018)`)
retail <- read_csv("inputs/skills/retail.csv") %>% mutate(NAICS = "44-45")
healthcare <- read_csv("inputs/skills/healthcare.csv") %>% mutate(NAICS = "62")

waiter.transfer <- read_csv("inputs/skills/waiter transfer.csv") %>% filter()

```

```{r}
qcew.preCOVID.LA <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22000")

## getting numbers by parish to combine for the netro
qcew.preCOVID.075 <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22075")
qcew.preCOVID.087 <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area ="22087")
qcew.preCOVID.089 <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22089")
qcew.preCOVID.093 <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22093")
qcew.preCOVID.095 <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22095")
qcew.preCOVID.103 <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22103")

## these parishes also get used on their own
qcew.preCOVID.OP <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22071")

qcew.preCOVID.JP <- blsQCEW(method = "Area",
        year = "2019",
        quarter = 3,
        area = "22051")

qcew.NOmsa.par <- bind_rows(qcew.preCOVID.JP,qcew.preCOVID.075,qcew.preCOVID.087,qcew.preCOVID.089,qcew.preCOVID.093,qcew.preCOVID.095,qcew.preCOVID.103,qcew.preCOVID.OP)

LA.ind.emp <- qcew.preCOVID.LA %>% 
  filter(agglvl_code == 54) %>% 
  group_by(industry_code) %>% 
  summarise(LAemp = sum(month3_emplvl, na.rm = TRUE))

NOmsa.ind.emp <- qcew.NOmsa.par %>% 
  filter(agglvl_code == 74) %>% 
  group_by(industry_code) %>% 
  summarise(NOmsaemp = sum(month3_emplvl, na.rm = TRUE))

```


### Expected number of unemployment claims by industry   

This table breaks down state-level initial + continued unemployment claims for the week ending June 6. The given industry's employment in the New Orleans metro is also shown, as is the expected number of UI claims in each industry. This number is calculated by multiplying the state-level rate of UI claims for the industry by the metro-level employment for the industry. The industries are arranged by their rank in terms of the percentnumber of employees expected to be claiming UI, and the rank by expected percent of metro UI claims is also shown.   

It's important to keep in mind that these extrapolated numbers are only helpful as a way to look at relative impacts between industries; they are not only underestimates, but also certainly inaccurate to some degree. See the final two graphs in our latest brief for more context on what this data actually means: https://www.datacenterresearch.org/reports_analysis/coming-into-focus-early-indicators-of-pandemic-job-loss-in-the-new-orleans-metro-area/


```{r}
lwc.ui.Industry.611 <- read_csv("inputs/LWC/lwcUIDashboard_20200611/Industry (Louisiana Unemployment Compensation Claims for Week Ending on March 28 2020 ).csv",
                            col_types = cols(`Week Ending` = col_date(format = "%m/%d/%Y"))) %>% 
  mutate(NAICS = str_replace(`NAICS Code`, "NAICS ","")) %>% 
  select(-`Number of Records`, -`NAICS Code`) %>% 
  mutate(Industry = ifelse(Industry == "Health Care and Social Assistence", "Health Care and Social Assistance", Industry),
         Industry = ifelse(Industry == "AccomModation & Food Services", "Accommodation & Food Services", Industry)) 

lwc.ui.Parish.611 <- read_csv("inputs/LWC/lwcUIDashboard_20200611/Parish (Louisiana Unemployment Compensation Claims for Week Ending on March 28 2020 ).csv", 
                          col_types = cols(`Week Ending` = col_date(format = "%m/%d/%Y"))) %>% 
    select(-`Number of Records`)

# NOmsa.claims <- lwc.ui.Parish.611 %>% 
#   filter(`Claim Type` == "Initial",
#          Parish %in% c("Orleans Parish",
#                                       "Jefferson Parish",
#                                       "St. Tammany Parish",
#                                       "St. John the Baptist Parish",
#                                       "St. Bernard Parish",
#                                       "Plaquemines Parish",
#                                       "St. Charles Parish",
#                                       "St. James Parish")) %>%
#   group_by(`Week Ending`) %>%
#   summarise(Claims = sum(Claims))

LA.ind.totclaims <- lwc.ui.Industry.611 %>% 
  group_by(NAICS, Industry, `Week Ending`) %>% 
  # filter(`Week Ending` >= "2020-03-21") %>% 
  summarise(LAUIclaims = sum(Claims))

NOmsa.ind.totclaims.ex <- LA.ind.totclaims %>% 
  left_join(LA.ind.emp, by=c("NAICS" = "industry_code")) %>% 
  filter(Industry != "Unclassified") %>% 
  mutate(LAUIclaims.p = LAUIclaims/LAemp) %>%
  left_join(NOmsa.ind.emp, by =c("NAICS" = "industry_code")) %>%
  mutate(expMSAclaims.ind = NOmsaemp * LAUIclaims.p) %>%
  filter(`Week Ending`==max(`Week Ending`)) %>%
  ungroup() %>%
  mutate(rank.p = dense_rank(-LAUIclaims.p),
         rank = dense_rank(-expMSAclaims.ind)) 

NOmsa.ind.totclaims.ex%>%
  select(Industry, pct_claiming_UI = LAUIclaims.p, emp_in_NOmsa = NOmsaemp, exp_num_claims_NOmsa = expMSAclaims.ind, rank_by_num_exp_claims = rank, rank_by_pct_claiming_UI = rank.p) %>%
  arrange(rank_by_num_exp_claims) %>%
  kableExtra::kable(digits = 2)
```

### Occupations with highest expected number of UI claims

After identifying the three industries losing the larget number of workers according to UI claims in the table above, I queried each industry's metro-level employment by occupation using "Staffing Patterns" in EMSI. I then combined the industries' occupation data, calculated the expected number of employees in each industry claiming UI as it was done for the chart above, and selected to top 30 occupations by number claiming UI.   

So, for example, at the state level, about 38% of people employed in the Acomodation & Food Services industry claimed unemployment in the week ending June 6. Within that industry is the occupation "Waiters and Waitresses," which employs about 14,300 people in the metro. If 38% of them claim UI, that's about 5,300 claims in that occupation.   


```{r}
occEmp <- NOmsa.ind.totclaims.ex %>%
  filter(NAICS %in% c("72","44-45", "62")) %>%
  right_join(bind_rows(food,retail,healthcare), by = "NAICS") %>%
  mutate(num.unemployed = LAUIclaims.p* `Employed in Industry (2018)`)%>% top_n(n = 30, wt = num.unemployed)

occEmp %>%
  select(Industry, Description, pct_claiming_UI = LAUIclaims.p, `Employed in Industry (2019)`, exp_num_claiming_UI = num.unemployed) %>%
  arrange(-exp_num_claiming_UI) %>%
  kable(digits = 2)
```


### Occupations that are most mathematically similar to "Waiters and Waitresses," based on O*NET skills scores.    

O*NET skills data describes each of 900+ occupations by the level & importance of a set of skills for carrying out the job, e.g. "Active Listening" and "Troubleshooting." The following table shows the handful of these 900+ occupations that are most similar to the occupation we estimate to have lost the highest number of jobs in the metro - waiters and waitresses. This similarity was identified using a simple kmeans clustering algorithm to divide the data into 50 clusters of similar occupations.   


```{r}
skills.wide <- skills %>%
  mutate(elID = paste0(`Element Name`,`Scale ID`)) %>%
  select(Title, elID, `Data Value`) %>%
  spread(key = elID, value = `Data Value`)

titles <- skills.wide %>% select(Title)

scores <- skills.wide %>% select(-Title)

set.seed(12)
k_skills <- kmeans((scores), centers = 50, nstart = 25) 

clust_skills <-skills.wide %>%
  mutate(cluster = k_skills$cluster)

wwc <- clust_skills %>%
  filter(Title == "Waiters and Waitresses") %>%
  pull(cluster)

clust_skills %>%
  filter(cluster == wwc) %>%
  kable()
```


### A note on technology skills

O*NET also provides data on the technology skills required for the 900+ occupations. This data is pretty different from the data set shown above in form -- it is just a list of every technology used in each occupation; e.g. "POS software" for waiters/waitresses and "Microsoft Excel" for CEOs. It would be easy to do a basic analysis of the technologies that people would need to be trained on in order to transfer from one occupation to another. 

