14 Validity and Causal Reasoning

Throughout this book, you have loaded data, summarized it, visualized it, and modeled relationships. This final content chapter addresses the most consequential question in data analysis: does the pattern you found reflect a real causal mechanism, or is it an artifact of how the data was collected, analyzed, or presented?

This chapter introduces types of research validity, counterfactual reasoning, and quasi-experimental designs. We will use a concrete example — whether cannabis legalization in Colorado reduced opioid overdose deaths — to illustrate how these concepts shape the way you analyze data and communicate findings. The analysis draws on data from the CDC’s WONDER database and the research design frameworks from Freese and Kevern’s Handbook of Causal Analysis for Social Research and Reichardt’s entry on internal validity in the International Encyclopedia of the Social & Behavioral Sciences.

14.1 Setup

import numpy as np

import pandas as pd
idx = pd.IndexSlice
pd.options.display.max_columns = 200

%matplotlib inline
import matplotlib.pyplot as plt
from matplotlib import ticker

import seaborn as sb
from scipy import stats

14.2 The Data

We use mortality data from the CDC’s WONDER database, which provides county-level death counts by cause from 1999 through 2020. We focus on drug/alcohol-induced causes using the “Unintentional” (X40–X44), “Suicide” (X60–X64), and “Undetermined” (Y10–Y14) ICD-10 codes to examine opioid overdose trends.

cod_df = pd.read_csv('wonder_cod.txt', sep='\t', dtype={'County Code': str})

# Extract the state from the "County" field (e.g., "Denver County, CO")
cod_df['State'] = cod_df['County'].str.split(', ').apply(lambda x: x[-1])

# Compute crude rate per 100k residents
cod_df['Rate'] = cod_df['Deaths'] / cod_df['Population'] * 100000

cod_df.head()

Not all counties are captured: the CDC suppresses small counts for privacy reasons, leaving only about 1,380 of the roughly 3,350 U.S. counties.

len(cod_df['County Code'].unique())

14.2.1 State-level rates

Pivot the data to compute annual state-level death rates per 100,000 residents.

deaths = pd.pivot_table(
    index='State', columns='Year', values='Deaths',
    aggfunc='sum', data=cod_df
)

population = pd.pivot_table(
    index='State', columns='Year', values='Population',
    aggfunc='sum', data=cod_df
)

state_rate_pivot = deaths / population * 100000
state_rate_pivot

14.3 The Question

Colorado legalized recreational cannabis on January 1, 2014. Some researchers have hypothesized that cannabis legalization reduces opioid overdose deaths because patients substitute cannabis for opioids in pain management. Let’s see what the data says.

legal_states = ['CO', 'WA']
nonlegal_states = [s for s in state_rate_pivot.index if s not in legal_states]

f, (ax0, ax1) = plt.subplots(1, 2, figsize=(12, 6))

ax0.set_ylabel('Deaths per 100k')
ax1.set_ylabel('Deaths per 100k')

state_rate_pivot.loc[legal_states, :].T.plot(ax=ax0, zorder=2, lw=2)
state_rate_pivot.loc[nonlegal_states].T.plot(ax=ax0, c='k', linewidth=1,
                                              alpha=.25, zorder=1, legend=False)

state_rate_pivot.loc[legal_states, :].T.plot(ax=ax1, zorder=2, lw=2)
state_rate_pivot.loc[nonlegal_states].mean(0).plot(ax=ax1, c='k',
                                                    linewidth=5, zorder=1, label='Non-legal avg')
ax1.legend()

f, ax = plt.subplots(1, 1, figsize=(12, 6))

state_rate_pivot.loc[legal_states, :].T.plot(ax=ax, zorder=2, lw=2)
state_rate_pivot.loc[nonlegal_states].mean(0).plot(ax=ax, c='k', linewidth=5, zorder=1)

ax.axvline(2014, color='r', ls='--', zorder=0)
f.tight_layout()

At first glance, Colorado’s rate appears to flatten around 2014 while the national average continues to climb. But is this evidence of a causal effect?

14.4 Types of Validity

14.4.1 Internal validity

Internal validity asks: did the treatment (cannabis legalization) cause the observed outcome (change in opioid deaths)? Threats to internal validity include confounding variables (other things that changed around 2014), selection bias (Colorado is not a random state), and maturation (trends that would have occurred anyway regardless of legalization).

14.4.2 External validity

External validity asks: would the same effect occur in other states? Even if legalization truly reduced opioid deaths in Colorado, the effect might not generalize to states with different demographics, healthcare systems, or drug markets.

14.4.3 Construct validity

Construct validity asks: are we measuring what we think we’re measuring? “Drug-induced deaths” is a broad category that includes many substances beyond opioids. Changes in coding practices, reporting requirements, or medical examiner staffing could affect the count without any real change in behavior.

14.5 Single-Group Pre/Post Design

The simplest design compares Colorado’s rates before and after legalization.

co_cod_df = cod_df[(cod_df['State'] == 'CO') & (cod_df['Year'] >= 2008)]

f, ax = plt.subplots(figsize=(12, 6))
sb.boxplot(x='Year', y='Rate', data=co_cod_df, ax=ax, color='green')
ax.axvline(5.5, color='r', ls='--', zorder=2, label='Legalization (2014)')
ax.axhline(
    co_cod_df.loc[co_cod_df['Year'] == 2014, 'Rate'].mean(),
    color='r', ls='--', zorder=2
)
ax.set_ylim((0, 50))
ax.legend()

co_before = co_cod_df.loc[co_cod_df['Year'].isin(range(2008, 2014)), 'Rate']
co_after = co_cod_df.loc[co_cod_df['Year'].isin(range(2015, 2021)), 'Rate']

print(f"Before mean: {co_before.mean():.1f}")
print(f"After mean: {co_after.mean():.1f}")
print(f"Difference: {co_after.mean() - co_before.mean():.1f}")

stats.ttest_ind(co_after, co_before, equal_var=False)

The fatal flaw of pre/post designs

Even if the difference is statistically significant, we cannot attribute it to legalization. Many other things changed between 2008–2013 and 2014–2020: the opioid crisis intensified nationally, fentanyl entered the drug supply, and Medicaid expansion occurred. Without a comparison group, we cannot distinguish the effect of legalization from the effect of time.

14.6 Multiple-Group Comparison Design

A stronger design compares Colorado to states that did not legalize cannabis but had similar pre-legalization death rates.

state_rate_pivot.loc['CO', 2008:2013].mean()

state_rate_pivot.loc[:, 2008:2013].mean(1).sort_values()

We identify comparison states with similar pre-2014 rates and track all three trajectories.

co_before_df = cod_df.loc[
    (cod_df['State'] == 'CO') & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate'
]
co_after_df = cod_df.loc[
    (cod_df['State'] == 'CO') & (cod_df['Year'] > 2014), 'Rate'
]

mt_before_df = cod_df.loc[
    (cod_df['State'] == 'MT') & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate'
]
mt_after_df = cod_df.loc[
    (cod_df['State'] == 'MT') & (cod_df['Year'] > 2014), 'Rate'
]

de_before_df = cod_df.loc[
    (cod_df['State'] == 'DE') & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate'
]
de_after_df = cod_df.loc[
    (cod_df['State'] == 'DE') & (cod_df['Year'] > 2014), 'Rate'
]

f, ax = plt.subplots(1, 1, figsize=(12, 6))

state_rate_pivot.loc[['CO', 'MT', 'DE', 'NM'], :].T.plot(ax=ax, zorder=2, lw=3)
state_rate_pivot.loc[nonlegal_states].mean(0).plot(ax=ax, c='k', linewidth=5, zorder=1)

ax.axvline(2014, color='r', ls='--', zorder=0)
f.tight_layout()

14.6.1 Visualize the distributions

f, axs = plt.subplots(
    nrows=1, ncols=3, figsize=(12, 4),
    subplot_kw={'xlim': (0, 100), 'ylim': (-.005, .2), 'xlabel': 'Opioid deaths per 100k'},
    sharey=True
)

co_before_df.plot.kde(ax=axs[0], color='blue')
co_after_df.plot.kde(ax=axs[0], color='red')
axs[0].axvline(co_before_df.mean(), color='blue', lw=.5)
axs[0].axvline(co_after_df.mean(), color='red', lw=.5)
axs[0].set_title('Colorado')

mt_before_df.plot.kde(ax=axs[1], color='blue')
mt_after_df.plot.kde(ax=axs[1], color='red')
axs[1].axvline(mt_before_df.mean(), color='blue', lw=.5)
axs[1].axvline(mt_after_df.mean(), color='red', lw=.5)
axs[1].set_title('Montana')

de_before_df.plot.kde(ax=axs[2], color='blue')
de_after_df.plot.kde(ax=axs[2], color='red')
axs[2].axvline(de_before_df.mean(), color='blue', lw=.5)
axs[2].axvline(de_after_df.mean(), color='red', lw=.5)
axs[2].set_title('Delaware')

f.tight_layout()

14.6.2 Compare changes across states

for name, before, after in [('Colorado', co_before_df, co_after_df),
                             ('Montana', mt_before_df, mt_after_df),
                             ('Delaware', de_before_df, de_after_df)]:
    diff = after.mean() - before.mean()
    t, p = stats.ttest_ind(after, before, equal_var=False)
    print(f"{name}: before {before.mean():.1f}, after {after.mean():.1f}, "
          f"diff {diff:+.1f}, t={t:.2f}, p={p:.4f}")

Montana’s increase is non-significant (p > 0.05) while Delaware’s increase is highly significant (p < 0.001), reflecting the fentanyl crisis that hit East Coast states particularly hard.

14.7 Counterfactual Reasoning

The core challenge of causal inference is that we can never observe the counterfactual — what would have happened to Colorado if it had not legalized cannabis. The comparison group serves as a proxy for this counterfactual, but it is only as good as the match between the treatment and comparison groups.

This is why randomized controlled experiments are the gold standard: random assignment ensures comparability on all characteristics, observed and unobserved. Observational studies like ours can never fully eliminate confounding.

14.8 Difference-in-Differences

A common quasi-experimental approach is difference-in-differences (DiD): compute the before-to-after change in the treatment group and subtract the before-to-after change in the comparison group.

co_change = co_after_df.mean() - co_before_df.mean()
mt_change = mt_after_df.mean() - mt_before_df.mean()

did = co_change - mt_change

print(f"Colorado change: {co_change:+.1f}")
print(f"Montana change: {mt_change:+.1f}")
print(f"Difference-in-differences: {did:+.1f}")

A negative DiD estimate would suggest that Colorado’s rate increased less (or decreased more) than Montana’s, consistent with a protective effect. But even this design relies on the parallel trends assumption — that Colorado and Montana would have followed the same trajectory in the absence of legalization.

14.9 Communicating Causal Claims

When you present findings from observational data, your language must reflect the limitations of your design.

Don’t write: “Cannabis legalization reduced opioid deaths in Colorado.”

Write instead: “Colorado’s drug-induced death rate increased less than comparable non-legalization states between 2008 and 2020, though this cross-state comparison cannot definitively attribute the difference to legalization.”

This is not hedging for its own sake — it is honest communication about what the evidence can and cannot support. As discussed in Chapter 11, your credibility depends on being clear about the boundaries of your knowledge.

Exercise 1: Identify threats to validity

List three specific threats to internal validity for the cannabis-opioid analysis. For each, explain how it could bias the results in either direction.

Solution:

Confounding with Medicaid expansion: Colorado expanded Medicaid in 2014, the same year as legalization. Increased access to addiction treatment could reduce overdose deaths independently of cannabis access, inflating the apparent protective effect.
Changes in drug supply: The rise of fentanyl accelerated after 2014, disproportionately affecting some states. If fentanyl reached Colorado later than comparison states, the apparent protective effect would be overstated.
Reporting changes: Improvements in medical examiner toxicology testing could change how deaths are classified. A state that began testing for more substances might reclassify deaths from “drug-induced” to more specific categories, artificially reducing the count.

Exercise 2: Choose better comparison states

Using the state_rate_pivot data, identify the three states whose pre-2014 death rates most closely match Colorado’s. Plot their trajectories alongside Colorado’s. Does the visual evidence support or undermine the legalization hypothesis?

Solution:

co_pre = state_rate_pivot.loc['CO', 2008:2013].mean()
distances = (state_rate_pivot.loc[:, 2008:2013].mean(1) - co_pre).abs()
best_matches = distances.drop(['CO', 'WA']).sort_values().head(3).index.tolist()

print(f"Best matches: {best_matches}")

f, ax = plt.subplots(figsize=(12, 6))
state_rate_pivot.loc[['CO'] + best_matches, :].T.plot(ax=ax, linewidth=2)
ax.axvline(2014, color='black', linestyle=':')
ax.set_ylabel('Deaths per 100,000')
ax.set_title('Colorado vs. Best-Matched Comparison States')
ax.legend()

Exercise 3: Compute a difference-in-differences

Using the comparison states you identified in Exercise 2, compute the difference-in-differences estimate for each pair (Colorado vs. each comparison state). Are the estimates consistent across comparisons?

Solution:

co_change = co_after_df.mean() - co_before_df.mean()

for state in best_matches:
    before = cod_df.loc[
        (cod_df['State'] == state) & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate'
    ]
    after = cod_df.loc[
        (cod_df['State'] == state) & (cod_df['Year'] > 2014), 'Rate'
    ]
    change = after.mean() - before.mean()
    did = co_change - change
    print(f"CO vs {state}: DiD = {did:+.1f}")

If the DiD estimates are consistent (all negative, or all positive, or all near zero), the finding is more robust. If they vary widely, the result is sensitive to the choice of comparison group — a red flag.

Exercise 4: Write a data brief

Write a 200-word data brief for a state legislator considering cannabis legalization. Include one visualization, one key statistic, and explicit acknowledgment of what the evidence can and cannot tell us.

Solution:

This is an open-ended exercise. A strong response leads with the finding, provides one clear chart, cites the DiD estimate, and concludes with caveats about confounding and the limits of cross-state comparison. The brief should avoid claiming causation while being clear about the direction and magnitude of the observed association.

14.10 Looking Back

This chapter concludes the content of the book. You have traveled from the mechanics of loading a CSV file to the philosophical foundations of causal inference. Along the way, you developed skills in data shaping, visualization, temporal analysis, statistical modeling, and — most importantly — the communication of findings in a way that is honest, clear, and appropriate for your audience.

These skills are not separate competencies to be checked off a list. They form an integrated practice of thinking carefully about data and the stories it can and cannot tell. Every chart you make embeds assumptions about what matters and what does not. Every regression coefficient you report implies a theory about how the world works. Every audience you address deserves an honest accounting of what you know and what remains uncertain.

The best data analysts are not the ones who produce the most sophisticated models or the most polished charts. They are the ones who ask the right questions, acknowledge the right limitations, and communicate with clarity and integrity.

# Validity and Causal Reasoning {#sec-validity} Throughout this book, you have loaded data, summarized it, visualized it, and modeled relationships. This final content chapter addresses the most consequential question in data analysis: **does the pattern you found reflect a real causal mechanism, or is it an artifact of how the data was collected, analyzed, or presented?** This chapter introduces types of research validity, counterfactual reasoning, and quasi-experimental designs. We will use a concrete example — whether cannabis legalization in Colorado reduced opioid overdose deaths — to illustrate how these concepts shape the way you analyze data and communicate findings. The analysis draws on data from the CDC's WONDER database and the research design frameworks from Freese and Kevern's *Handbook of Causal Analysis for Social Research* and Reichardt's entry on internal validity in the *International Encyclopedia of the Social & Behavioral Sciences*. ## Setup ```{python} import numpy as np import pandas as pd idx = pd.IndexSlice pd.options.display.max_columns = 200 %matplotlib inline import matplotlib.pyplot as plt from matplotlib import ticker import seaborn as sb from scipy import stats ``` ## The Data We use mortality data from the CDC's [WONDER](https://wonder.cdc.gov/) database, which provides county-level death counts by cause from 1999 through 2020. We focus on drug/alcohol-induced causes using the "Unintentional" (X40–X44), "Suicide" (X60–X64), and "Undetermined" (Y10–Y14) ICD-10 codes to examine opioid overdose trends. ```{python} cod_df = pd.read_csv('wonder_cod.txt', sep='\t', dtype={'County Code': str}) # Extract the state from the "County" field (e.g., "Denver County, CO") cod_df['State'] = cod_df['County'].str.split(', ').apply(lambda x: x[-1]) # Compute crude rate per 100k residents cod_df['Rate'] = cod_df['Deaths'] / cod_df['Population'] * 100000 cod_df.head() ``` Not all counties are captured: the CDC suppresses small counts for privacy reasons, leaving only about 1,380 of the roughly 3,350 U.S. counties. ```{python} len(cod_df['County Code'].unique()) ``` ### State-level rates Pivot the data to compute annual state-level death rates per 100,000 residents. ```{python} deaths = pd.pivot_table( index='State', columns='Year', values='Deaths', aggfunc='sum', data=cod_df ) population = pd.pivot_table( index='State', columns='Year', values='Population', aggfunc='sum', data=cod_df ) state_rate_pivot = deaths / population * 100000 state_rate_pivot ``` ## The Question Colorado legalized recreational cannabis on January 1, 2014. Some researchers have hypothesized that cannabis legalization reduces opioid overdose deaths because patients substitute cannabis for opioids in pain management. Let's see what the data says. ```{python} legal_states = ['CO', 'WA'] nonlegal_states = [s for s in state_rate_pivot.index if s not in legal_states] ``` ```{python} f, (ax0, ax1) = plt.subplots(1, 2, figsize=(12, 6)) ax0.set_ylabel('Deaths per 100k') ax1.set_ylabel('Deaths per 100k') state_rate_pivot.loc[legal_states, :].T.plot(ax=ax0, zorder=2, lw=2) state_rate_pivot.loc[nonlegal_states].T.plot(ax=ax0, c='k', linewidth=1, alpha=.25, zorder=1, legend=False) state_rate_pivot.loc[legal_states, :].T.plot(ax=ax1, zorder=2, lw=2) state_rate_pivot.loc[nonlegal_states].mean(0).plot(ax=ax1, c='k', linewidth=5, zorder=1, label='Non-legal avg') ax1.legend() ``` ```{python} f, ax = plt.subplots(1, 1, figsize=(12, 6)) state_rate_pivot.loc[legal_states, :].T.plot(ax=ax, zorder=2, lw=2) state_rate_pivot.loc[nonlegal_states].mean(0).plot(ax=ax, c='k', linewidth=5, zorder=1) ax.axvline(2014, color='r', ls='--', zorder=0) f.tight_layout() ``` At first glance, Colorado's rate appears to flatten around 2014 while the national average continues to climb. But is this evidence of a causal effect? ## Types of Validity ### Internal validity Internal validity asks: **did the treatment (cannabis legalization) cause the observed outcome (change in opioid deaths)?** Threats to internal validity include confounding variables (other things that changed around 2014), selection bias (Colorado is not a random state), and maturation (trends that would have occurred anyway regardless of legalization). ### External validity External validity asks: **would the same effect occur in other states?** Even if legalization truly reduced opioid deaths in Colorado, the effect might not generalize to states with different demographics, healthcare systems, or drug markets. ### Construct validity Construct validity asks: **are we measuring what we think we're measuring?** "Drug-induced deaths" is a broad category that includes many substances beyond opioids. Changes in coding practices, reporting requirements, or medical examiner staffing could affect the count without any real change in behavior. ## Single-Group Pre/Post Design The simplest design compares Colorado's rates before and after legalization. ```{python} co_cod_df = cod_df[(cod_df['State'] == 'CO') & (cod_df['Year'] >= 2008)] f, ax = plt.subplots(figsize=(12, 6)) sb.boxplot(x='Year', y='Rate', data=co_cod_df, ax=ax, color='green') ax.axvline(5.5, color='r', ls='--', zorder=2, label='Legalization (2014)') ax.axhline( co_cod_df.loc[co_cod_df['Year'] == 2014, 'Rate'].mean(), color='r', ls='--', zorder=2 ) ax.set_ylim((0, 50)) ax.legend() ``` ```{python} co_before = co_cod_df.loc[co_cod_df['Year'].isin(range(2008, 2014)), 'Rate'] co_after = co_cod_df.loc[co_cod_df['Year'].isin(range(2015, 2021)), 'Rate'] print(f"Before mean: {co_before.mean():.1f}") print(f"After mean: {co_after.mean():.1f}") print(f"Difference: {co_after.mean() - co_before.mean():.1f}") stats.ttest_ind(co_after, co_before, equal_var=False) ``` ::: {.callout-warning} ## The fatal flaw of pre/post designs Even if the difference is statistically significant, we cannot attribute it to legalization. Many other things changed between 2008–2013 and 2014–2020: the opioid crisis intensified nationally, fentanyl entered the drug supply, and Medicaid expansion occurred. Without a comparison group, we cannot distinguish the effect of legalization from the effect of time. ::: ## Multiple-Group Comparison Design A stronger design compares Colorado to states that did not legalize cannabis but had similar pre-legalization death rates. ```{python} state_rate_pivot.loc['CO', 2008:2013].mean() ``` ```{python} state_rate_pivot.loc[:, 2008:2013].mean(1).sort_values() ``` We identify comparison states with similar pre-2014 rates and track all three trajectories. ```{python} co_before_df = cod_df.loc[ (cod_df['State'] == 'CO') & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate' ] co_after_df = cod_df.loc[ (cod_df['State'] == 'CO') & (cod_df['Year'] > 2014), 'Rate' ] mt_before_df = cod_df.loc[ (cod_df['State'] == 'MT') & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate' ] mt_after_df = cod_df.loc[ (cod_df['State'] == 'MT') & (cod_df['Year'] > 2014), 'Rate' ] de_before_df = cod_df.loc[ (cod_df['State'] == 'DE') & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate' ] de_after_df = cod_df.loc[ (cod_df['State'] == 'DE') & (cod_df['Year'] > 2014), 'Rate' ] ``` ```{python} f, ax = plt.subplots(1, 1, figsize=(12, 6)) state_rate_pivot.loc[['CO', 'MT', 'DE', 'NM'], :].T.plot(ax=ax, zorder=2, lw=3) state_rate_pivot.loc[nonlegal_states].mean(0).plot(ax=ax, c='k', linewidth=5, zorder=1) ax.axvline(2014, color='r', ls='--', zorder=0) f.tight_layout() ``` ### Visualize the distributions ```{python} f, axs = plt.subplots( nrows=1, ncols=3, figsize=(12, 4), subplot_kw={'xlim': (0, 100), 'ylim': (-.005, .2), 'xlabel': 'Opioid deaths per 100k'}, sharey=True ) co_before_df.plot.kde(ax=axs[0], color='blue') co_after_df.plot.kde(ax=axs[0], color='red') axs[0].axvline(co_before_df.mean(), color='blue', lw=.5) axs[0].axvline(co_after_df.mean(), color='red', lw=.5) axs[0].set_title('Colorado') mt_before_df.plot.kde(ax=axs[1], color='blue') mt_after_df.plot.kde(ax=axs[1], color='red') axs[1].axvline(mt_before_df.mean(), color='blue', lw=.5) axs[1].axvline(mt_after_df.mean(), color='red', lw=.5) axs[1].set_title('Montana') de_before_df.plot.kde(ax=axs[2], color='blue') de_after_df.plot.kde(ax=axs[2], color='red') axs[2].axvline(de_before_df.mean(), color='blue', lw=.5) axs[2].axvline(de_after_df.mean(), color='red', lw=.5) axs[2].set_title('Delaware') f.tight_layout() ``` ### Compare changes across states ```{python} for name, before, after in [('Colorado', co_before_df, co_after_df), ('Montana', mt_before_df, mt_after_df), ('Delaware', de_before_df, de_after_df)]: diff = after.mean() - before.mean() t, p = stats.ttest_ind(after, before, equal_var=False) print(f"{name}: before {before.mean():.1f}, after {after.mean():.1f}, " f"diff {diff:+.1f}, t={t:.2f}, p={p:.4f}") ``` Montana's increase is non-significant (p > 0.05) while Delaware's increase is highly significant (p < 0.001), reflecting the fentanyl crisis that hit East Coast states particularly hard. ## Counterfactual Reasoning The core challenge of causal inference is that we can never observe the **counterfactual** — what *would have happened* to Colorado if it had *not* legalized cannabis. The comparison group serves as a proxy for this counterfactual, but it is only as good as the match between the treatment and comparison groups. This is why randomized controlled experiments are the gold standard: random assignment ensures comparability on all characteristics, observed and unobserved. Observational studies like ours can never fully eliminate confounding. ## Difference-in-Differences A common quasi-experimental approach is **difference-in-differences** (DiD): compute the before-to-after change in the treatment group and subtract the before-to-after change in the comparison group. ```{python} co_change = co_after_df.mean() - co_before_df.mean() mt_change = mt_after_df.mean() - mt_before_df.mean() did = co_change - mt_change print(f"Colorado change: {co_change:+.1f}") print(f"Montana change: {mt_change:+.1f}") print(f"Difference-in-differences: {did:+.1f}") ``` A negative DiD estimate would suggest that Colorado's rate increased less (or decreased more) than Montana's, consistent with a protective effect. But even this design relies on the **parallel trends assumption** — that Colorado and Montana would have followed the same trajectory in the absence of legalization. ## Communicating Causal Claims When you present findings from observational data, your language must reflect the limitations of your design. **Don't write:** "Cannabis legalization reduced opioid deaths in Colorado." **Write instead:** "Colorado's drug-induced death rate increased less than comparable non-legalization states between 2008 and 2020, though this cross-state comparison cannot definitively attribute the difference to legalization." This is not hedging for its own sake — it is honest communication about what the evidence can and cannot support. As discussed in @sec-uncertainty, your credibility depends on being clear about the boundaries of your knowledge. ::: {.callout-note collapse="true"} ## Exercise 1: Identify threats to validity List three specific threats to internal validity for the cannabis-opioid analysis. For each, explain how it could bias the results in either direction. **Solution:** 1. **Confounding with Medicaid expansion:** Colorado expanded Medicaid in 2014, the same year as legalization. Increased access to addiction treatment could reduce overdose deaths independently of cannabis access, inflating the apparent protective effect. 2. **Changes in drug supply:** The rise of fentanyl accelerated after 2014, disproportionately affecting some states. If fentanyl reached Colorado later than comparison states, the apparent protective effect would be overstated. 3. **Reporting changes:** Improvements in medical examiner toxicology testing could change how deaths are classified. A state that began testing for more substances might reclassify deaths from "drug-induced" to more specific categories, artificially reducing the count. ::: ::: {.callout-note collapse="true"} ## Exercise 2: Choose better comparison states Using the `state_rate_pivot` data, identify the three states whose pre-2014 death rates most closely match Colorado's. Plot their trajectories alongside Colorado's. Does the visual evidence support or undermine the legalization hypothesis? **Solution:** ```{python} co_pre = state_rate_pivot.loc['CO', 2008:2013].mean() distances = (state_rate_pivot.loc[:, 2008:2013].mean(1) - co_pre).abs() best_matches = distances.drop(['CO', 'WA']).sort_values().head(3).index.tolist() print(f"Best matches: {best_matches}") f, ax = plt.subplots(figsize=(12, 6)) state_rate_pivot.loc[['CO'] + best_matches, :].T.plot(ax=ax, linewidth=2) ax.axvline(2014, color='black', linestyle=':') ax.set_ylabel('Deaths per 100,000') ax.set_title('Colorado vs. Best-Matched Comparison States') ax.legend() ``` ::: ::: {.callout-note collapse="true"} ## Exercise 3: Compute a difference-in-differences Using the comparison states you identified in Exercise 2, compute the difference-in-differences estimate for each pair (Colorado vs. each comparison state). Are the estimates consistent across comparisons? **Solution:** ```{python} co_change = co_after_df.mean() - co_before_df.mean() for state in best_matches: before = cod_df.loc[ (cod_df['State'] == state) & (cod_df['Year'] >= 2008) & (cod_df['Year'] < 2014), 'Rate' ] after = cod_df.loc[ (cod_df['State'] == state) & (cod_df['Year'] > 2014), 'Rate' ] change = after.mean() - before.mean() did = co_change - change print(f"CO vs {state}: DiD = {did:+.1f}") ``` If the DiD estimates are consistent (all negative, or all positive, or all near zero), the finding is more robust. If they vary widely, the result is sensitive to the choice of comparison group — a red flag. ::: ::: {.callout-note collapse="true"} ## Exercise 4: Write a data brief Write a 200-word data brief for a state legislator considering cannabis legalization. Include one visualization, one key statistic, and explicit acknowledgment of what the evidence can and cannot tell us. **Solution:** This is an open-ended exercise. A strong response leads with the finding, provides one clear chart, cites the DiD estimate, and concludes with caveats about confounding and the limits of cross-state comparison. The brief should avoid claiming causation while being clear about the direction and magnitude of the observed association. ::: ## Looking Back This chapter concludes the content of the book. You have traveled from the mechanics of loading a CSV file to the philosophical foundations of causal inference. Along the way, you developed skills in data shaping, visualization, temporal analysis, statistical modeling, and — most importantly — the communication of findings in a way that is honest, clear, and appropriate for your audience. These skills are not separate competencies to be checked off a list. They form an integrated practice of thinking carefully about data and the stories it can and cannot tell. Every chart you make embeds assumptions about what matters and what does not. Every regression coefficient you report implies a theory about how the world works. Every audience you address deserves an honest accounting of what you know and what remains uncertain. The best data analysts are not the ones who produce the most sophisticated models or the most polished charts. They are the ones who ask the right questions, acknowledge the right limitations, and communicate with clarity and integrity.