```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
setwd("~/Github/bigdata/lectures/l15/")
library(tidyverse)
```

## Today's Class {.smaller}

1. Review
    - Winsorizing
    - Bayes
    - Predictions Contest 2
2. Observational Data
3. RCTs: Average Treatment Effects
4. RCTs: Heterogenous Treatment Effects
5. RCTs: Targeting?
6. HW5 Review?

# Review

## Winsorizing

The standard routine does the following. 

1. Pick some quantile (e.g. 1%). 
2. Find that quantile -- e.g. find the 99%ile abs(logerror)
3. Set values higher than that quantile to that quantile. 

## Bayes Rule

$$ P[\beta | X] \propto P[X|\beta]P[\beta] $$

- $P[\beta|X]$ is referred to as the posterior
- $P[X|\beta]$ is the *likelihood* (we've seen before)
- $P[\beta]$ is the prior

**Systematizing Uncertainty**

## Primary Sources of Uncertainty in Modeling

1. Within the model, there is uncertainty. 
    - What are parameter values
    - What value will an observation take
2. We are not certain if the model is correct
    - Should we use a different model
    - Should we average with a different model
3. We don't know if the data is correct
    - Is there a clerical error?
    - Is there measurement error?
        - E.g. Covid Case counts in the US last March -- woefully undercounted.

## Predictions 2 - Daily Case Counts

![](../../predictions/preds2-p4.png)

## Predictions 2 - Sundays only

![](../../predictions/preds2-p5.png)

## Predictions 2 - Actual Predictions

![](../../predictions/preds2-p1.png)

## Predictions 2 - MSE Winners

![](../../predictions/preds2-p2.png)

## Predictions 2 - Price-is-Right Winners

![](../../predictions/preds2-p3.png)

## Predictions 2 - Wrapup

- Actual result: **27655** [(link)](https://codowd.com/bigdata/predictions/pred2-result.png)
- MSE Winners: Caden Kalinowski, Cagdas Okay, Matias Pietruszka
- PIR Winners: Frank Li, Caden, Amy Maldonado

![](../../predictions/pred2_reg.png)

## Predictions 3

New competition posted (this morning).

- Purely optional. Not graded. No canvas component.

US Sunday case counts _this sunday_. Two numbers:

1. Prediction
2. P[20% prediction error]

Possibly 1 more _optional_ competition. Would likely run into the summer. 

# Observational Data

---

We have a question:

- Do job training programs for the unemployed improve worker's earnings?

We have data: 

- Job training uptake, future earnings

Can we answer the question?

## Formalizing This

<!-- flowcharts made using at https://flowchart.fun -->

![](dag1.png)

## Formalizing This

![](dag1b.png)

## Formalizing This

![](dag1a.png)


## In Context

![](dag1context.png)

- Either Path A or Path B can't exist, for every other possible variable U. 

## Problem

If we have either (possibly unobserved) variables *U* that affect both X and Y, or Y might affect X, identifying causal effects of X on Y is _impossible_ with your data. 

Why? We can't disentangle the separate possible sources of changes.

E.g. Suppose extremely hard working individuals are more likely to _apply_ for job training, and have higher future incomes. Or any other of a number of possible suppositions.

## Assumption

In observational data, it is hard to guarantee "*No other variables have causal effects on both X and Y*" and that "*Y doesn't drive X*".

- We will see some settings, on Thursday, where we can make these assumptions more reasonable.
    - IV
    - RD
    - Diff-in-Diff
    
- But they will still be _assumptions_.
    
## RCTs

For now, we turn to randomized controlled trials. 

RCTs solve this problem, by introducing a variable which selects *X*, without having any relationship with *Y* -- a random number. 

## RCTs

![](dag2.png)

## RCTs

![](dag2context.png)
    
## RCTs

RCTs work because they force a treatment that is _unrelated_ to any other possible variable, and which isn't being caused by *Y*. 

So we can rule out the things we were worried about.


## Job Training RCT

In the early 90s, a federal job training program ("JTPA") ran an RCT where some people (~20k) were given offers immediately and some had eligibility delayed by 18 months _at random_. 

There are a wide variety of outputs you could look at. We will look at the effect on income over the following 30 months on adults (~11k people). 

First, some summary stats

## JTPA Summary Stats

```{r,echo=F}
load("jtpa/jtpa.RData")
jtpa = jtpa %>% select(-Index)
```
```{r, message=F}
jtpa %>% select(-c(9:14)) %>% 
    group_by(offer) %>% 
    summarize(across(everything(),mean))
```

## Summary Stats 2

```{r, message=F}
jtpa %>% select(c(2,9:14)) %>% 
    group_by(offer) %>% 
    summarize(across(everything(),mean))
```

## Comments

There are some weird things going on. 
```{r}
summary(as.factor(jtpa$age))
summary(as.factor(jtpa$married))
summary(as.factor(jtpa$wkless13))
```


## Comments 2

"Treatment" here is an offer for access to the JTPA.

- Uptake (train) is not identical to offer of treatment.

```{r, message=F}
jtpa %>% group_by(offer) %>%
    summarize(train = mean(train))
```

## Comments 3

- *What treatment* **IS** is not uniform across genders

```{r, message=F}
jtpa %>% group_by(male) %>%
    select(male,OJT_JSA,classroom,train) %>%
    summarize(across(everything(),mean))
```

- OJT: on-job-training, JSA: job-search assistance, classroom - skills training. 

## Observational Data -- Mean Difference *on uptake*

- We could look at difference between those who use and don't use program. 
    - Remember -- there may be selection driving this.
    
```{r}
smod = summary(lm(y~train,data=jtpa))
signif(smod$coefficients,3)
```

## RCTs -- Mean Difference

- Take the difference in means between group assigned treatment and group assigned control.
    - Average Treatment Effect: effect of _offer_ (AKA intention-to-treat) 
    - We can test for difference from 0, build CIs, etc - for this difference in means.
    
```{r}
smod = summary(lm(y~offer,data=jtpa))
signif(smod$coefficients,3)
```   

## RCT - complications

Basic Mean difference is unbiased on average across universes where you ran this experiment.

- But what if there was some residual variation in other variables. 
    - Like, people who were AFDC recipients were more likely to recieve treatment, *by chance*.
- Alternately, what if we are interested in the effect on some subpopulation? E.g. Gender differences?
- Finally, uptake in treatment group was like 60%. So the treatment effect must be larger for people _who used the program_. Can we figure that out?


## RCT  -- Subgroup Analysis

```{r}
smod = summary(lm(y~offer*male,data=jtpa))
signif(smod$coefficients,3)
```


## RCT -- Basic Controls

```{r}
smod = summary(lm(y~.,data=jtpa[,-c(3,12,11,13)]))
signif(smod$coefficients,3)
```

- Interactions? What is the model we should choose?

## RCTs - ATT

- We could try to look at the average treatment effect on the treated. 

```{r}
uptake = mean(jtpa$train[jtpa$offer == 1]) # ~64%
smod = summary(lm(y~offer,data=jtpa))$coefficients
smod[2,1]/uptake #Scale up Coef by uptake rate. 
```

- Uncertainty there is not just uncertainty internal to linear model -- but also around uptake. Bootstrap for CIs if you like. 

## RCTs - Big Data?

How is this related to this course?

1. Entirely about predictions. 
    - Obscured by simplicity of analysis
2. That was *Average* treatment effects. What about individual TEs?
    - But effects are heterogenous (we saw that in gender)
    - We can estimate the TE for any given individual. (Some may be negative)
    - This gives us Conditional Average Treatment Effects (CATEs)
3. We may want to build *targeting policies* based on RCTs. 
    - E.g. I run an RCT of some ad campaign. 
        1. Run RCT
        2. Identify CATEs
        3. Find _optimal_ policy for targeting ads. 
        4. Profit
    

## TEs are predictions? 

$$\widehat{ATE} = \left[ \frac{1}{n_1} \sum_{i:~T=1} y_i\right] - \left[\frac{1}{n_0} \sum_{i: ~T=0} y_i\right] = \bar y_1- \bar y_0$$



But suppose we define each individual's treatment effect as the sum of their observed outcome, and the *unobserved* counterfactual outcome where their treatment status was flipped:

$$ \widehat{TE}_i = y_i(1)-y_i(0)$$
To do this, we would need to make a prediction about $y_i$ in the unobserved counterfactual. 

## Prediction in Counterfactual

Suppose for the counterfactual, we predict the mean outcome from the other group. 

$$\hat{y}_i(T_i-1) = \bar y_{T_i-1}$$
Then (for a treated observation $i$), we have:

$$ \widehat{TE}_i = y_i - \bar y_0$$
And for an untreated observation, similarly:

$$\widehat{TE}_i = \bar y_1 - y_i$$

## Similarity

We could then, form a new estimate of the ATE:

$$\widetilde{ATE} = \frac1n \sum_{i=1}^n \widehat{TE}_i$$

With a lot of rewriting of sums -- we can show that:

$$ \widetilde{ATE} = \widehat{ATE}$$

## Proof in Data

```{r}
ybar0 = mean(jtpa$y[jtpa$offer == 0])
ybar1 = mean(jtpa$y[jtpa$offer == 1])
ybar1-ybar0
jtpa_est = jtpa %>% 
    mutate(y1 = y*offer+(1-offer)*ybar1,
           y0 = y*(1-offer)+offer*ybar0) %>%
    mutate(TE = y1-y0)
mean(jtpa_est$TE)
```


## TEs as predictions: A brief diversion

How does this affect our interpretation of 'controls' and of 'subgroup analysis'? 

- We can do the same basic exercise, where 'controls' or subgroup analyses affect _our counterfactual predictions_ and rewrite our ATE estimate as a mean of individual treatment effects. 

## TEs as predictions

First off, purely for calculating ATEs, this suggests a simple heuristic:

- If we can improve our OOS predictions for either treatment or control, we can improve our ATE estimate
    - We've seen a lot of ways one could improve a prediction in this course. 
    
## ATE from predictions:
```{r,echo=F}
library(ranger)
jtpa_full = jtpa
jtpa = jtpa[,-c(3,11:13)] #dropping "train" and other details of the treatment. Probably useful, but interpretation is weird.
```
```{r,message=F}
#Build Treat and control dfs
jtpa_cont = jtpa %>% filter(offer == 0)
jtpa_treat = jtpa %>% filter(offer == 1) 
#Estimate treat and control models
mod_treat = ranger(y~.-offer,data=jtpa_treat)
mod_cont = ranger(y~.-offer,data=jtpa_cont)
#Predict counterfactuals for data from other model
jtpa_cont$confact = predict(mod_treat,data = jtpa_cont)$predictions
jtpa_treat$confact  = predict(mod_cont,data = jtpa_treat)$predictions
#Estimate TEs
jtpa_cont$TE = jtpa_cont$confact - jtpa_cont$y
jtpa_treat$TE = jtpa_treat$y - jtpa_treat$confact
#Recombine
jtpa_est = rbind(jtpa_cont,jtpa_treat)
mean(jtpa_est$TE) #ATE
```

## ATE from Predictions

We could use _any_ of the many routines we've seen, and the many more you will encounter for our _counterfactual predictions_. 

- KNN
- Forests
- Logit
- LASSO
- Boosted Models
- Bagged Models
- Other Ensemble Models

Once this is a prediction problem -- we can do a lot. 

## TEs

But we could also use those TE estimates from individuals. 

```{r, echo=F}
jtpa_est = jtpa_est %>% mutate(y1 = y*offer+(1-offer)*confact,
                               y0 = y*(1-offer)+offer*confact)
ggplot(jtpa_est[sample(nrow(jtpa_est)),],
       aes(x=y0,y=y1,col=offer == 1)) +
    geom_point(alpha=0.05) + 
    geom_abline(intercept=0,slope=1,lty=2) 
```

## Individual TEs

```{r,echo=F}
ggplot(jtpa_est[sample(nrow(jtpa_est)),],
       aes(x=y0,y=y1)) +
    geom_hex()
```

## Individual TEs

```{r,echo=F}
ggplot(jtpa_est,aes(x=TE)) +
    geom_histogram(bins=100) +
    geom_vline(xintercept = mean(jtpa_est$TE),col=2,lty=2)
```

There is a lot of variation around our mean of `r mean(jtpa_est$TE)` in benefits. The standard deviation is `r sd(jtpa_est$TE)`. 

---

Is this variation in individual TEs sensible? 

- Probably? \$50k effects either way seem large, but only a small portion, `r mean(abs(jtpa_est$TE)>50000)`, experience that scale.
- Intuitively, some individuals got training that was beneficial, and some spent time doing training that they could have spent _working_ at other jobs that pay. So some loss seems about right. 

## Individual TEs, Do we care?

It looks like some fraction of individuals lost ~\$50k by engaging in this program. Not to mention the program cost to the government. 

- What if we could target the program to people who benefit?
- In other settings, like marketing, we may wish to target groups for whom the expense of advertising is less than the gain in revenue from those individuals.

$\implies$ Targeting. Can we use the RCT data for targeting?

## Targeting

We estimated treatment effects for individuals in an experiment.

- This means we had to predict the unobserved counterfactual
- But it also means we witnessed the outcome for one treatment possibility.

To do targeting we need to:

1) Make predictions for individuals under treatment
2) Make predictions for individuals under control
3) Estimate TEs for each individual
4) Compare TEs to some threshold (\$0? \$1k?) to determine eligibility.

## Targeting

1. Fully a counterfactual exercise.
2. We need _good_ OOS predictions
    - We may want to do some CV to determine performance
3. We need something to compare to. "Opportunity Costs" 

Talk more Thursday. 

## An Aside: RCTs and Mechanisms

RCTs are not good at identifying the mechanisms for effects. 

This is because identifying a specific mechanism means establishing a link from X to that mechanism, from that mechanism to Y, **and** the lack of links between X and Y through other mechanisms.

You need multiple overlapping RCTs to do this kind of thing -- each looking at different things. 

## An Aside: RCTs and Mechanisms

![](dag_mechanism.png)



    
# HW 5 Review (if time)

# Wrap up

## Things to do

Homework 6 is due tomorrow. New prediction competition was released yesterday -- purely optional. 

On Thursday we will:

1. Wrap up targeting/RCTs
2. *Briefly* encounter a few 'observational' causal methods:
    a. Instrumental Variables
        a. Useful for "Intention to Treat" vs "Average effect on Treated"
    b. Regression Discontinuity
    c. Diff-in-Diff
        a. SCM

# Bye!








