Empirical Research Methods

Research Code Showcase

Sample Stata code from a study examining food insecurity before and after gray divorce. The workflow demonstrates panel-data preparation, Difference-in-Differences estimation, propensity score matching, marital-transition analysis, and event-study visualization.

Stata Causal Inference Updated July 2026
Stata research code and statistical analysis

Project Overview

This demonstration presents a simplified empirical workflow for analyzing the relationship between gray divorce and food insecurity using longitudinal survey data. The code illustrates common methods used in applied economics, health economics, and public-policy research.

  • Panel Data
  • Data Reshaping
  • Difference-in-Differences
  • Fixed Effects
  • Propensity Score Matching
  • Event Study
  • Clustered Standard Errors

Section 01

Panel-Data Preparation

The first step loads the Health and Retirement Study panel, retains the variables required for analysis, reshapes the data from wide to long format, assigns survey years, and constructs marital-status and income variables.

Data loading, reshaping, and variable construction Stata 
* Gray Divorce Study - Food Insecurity Analysis
* Load and reshape HRS panel data

use "HRS_grayDiv.dta", clear

keep hhidpn foodInsec* female ragey_b rmstat ///
    hhIncome baseline_weight

reshape long foodInsec@ ragey_b@ rmstat@ hhIncome@, ///
    i(hhidpn) j(period) string


* Create survey-year variable

gen year = .

replace year = 1998 if period == "4"
replace year = 2000 if period == "5"
replace year = 2002 if period == "6"


* Define marital status

gen marital_s = .

replace marital_s = 1 if rmstat >= 1 & rmstat <= 3
replace marital_s = 2 if rmstat >= 4 & rmstat <= 6

label define marital_s ///
    1 "Married" ///
    2 "Divorced"

label values marital_s marital_s


* Create treatment and post-divorce indicators

bysort hhidpn (year): ///
    egen everDiv = max(marital_s == 2)

gen treat = everDiv
gen post  = marital_s == 2


* Transform household income

gen loghhincome = log(hhIncome)
gen age_sq      = ragey_b^2
Purpose: Reshaping the data creates one observation per individual and survey wave, making the dataset suitable for longitudinal fixed-effects and event-study models.

Section 02

Difference-in-Differences Analysis

The model estimates the change in food insecurity associated with divorce while controlling for age, household income, survey-year effects, and time-invariant individual characteristics. Standard errors are clustered at the individual level.

Individual fixed-effects DiD model Stata 
* Difference-in-Differences regression

areg foodInsec ///
    i.treat##i.post ///
    ragey_b ///
    age_sq ///
    loghhincome ///
    i.year ///
    [pweight = baseline_weight], ///
    absorb(hhidpn) ///
    vce(cluster hhidpn)
Interpretation: The coefficient on the treatment-by-post interaction measures the estimated change in food insecurity among individuals experiencing divorce, relative to the comparison group.

Section 03

Propensity Score Matching

Propensity scores are estimated using observable demographic and socioeconomic characteristics. The score distributions are examined before one-to-one nearest-neighbor matching is conducted with a caliper and without replacement.

Propensity-score estimation and diagnostics Stata 
* Step 1: Estimate propensity scores

local treatment "treat"

local match_var ///
    "female ragey_b loghhincome hispanic black othrace
     leshsdg hsdeg college snap"

logit `treatment' `match_var' ///
    if `treatment' != . & n == 1

predict ps if n == 1, pr


* Step 2: Summarize propensity-score distributions

summarize ps if n == 1

tabstat ps, ///
    statistics(mean sd) ///
    by(`treatment')


* Plot overlapping propensity-score distributions

twoway ///
    (histogram ps if `treatment' == 0 & n == 1, ///
        bin(50) start(0) ///
        lcolor(gray) fcolor(gray)) ///
    (histogram ps if `treatment' == 1 & n == 1, ///
        bin(50) start(0) ///
        lcolor(black) fcolor(none)), ///
    legend( ///
        order(1 "Stable Partnership" ///
              2 "Approaching Divorce") ///
        position(12) cols(1)) ///
    graphregion( ///
        fcolor(white) ///
        ifcolor(white) ///
        ilcolor(white)) ///
    xtitle("Propensity Score") ///
    xlabel(0(0.1)0.4, format(%3.1f)) ///
    ytitle("Density of Observations") ///
    xscale(noline) ///
    yscale(noline)

One-to-one matching and matched-pair construction Stata 
* Step 3: One-to-one nearest-neighbor matching

psmatch2 treat if n == 1, ///
    pscore(ps) ///
    outcome(foodInsec) ///
    caliper(0.001) ///
    noreplacement ///
    neighbor(1)


* Construct matched-pair identifiers

gen pair = _id if _treated == 0 & n == 1

replace pair = _n1 ///
    if _treated == 1 & n == 1

bysort pair: ///
    egen paircount = count(pair) if n == 1

gen validpair = ///
    paircount > 1 & n == 1

bysort hhidpn: ///
    egen validpair_all = max(validpair)

drop if validpair_all != 1

drop pair paircount validpair validpair_all
Research consideration: Matching quality should also be evaluated through covariate-balance statistics and common-support diagnostics, rather than relying only on propensity-score overlap.

Section 04

Remarriage and Marital Transitions

This section identifies changes in marital status over time, counts marital transitions, and constructs an indicator for individuals who experience remarriage after divorce.

Longitudinal marital-transition indicators Stata 
* Declare the panel structure

xtset hhidpn year

gen marital_r = marital_s


* Identify changes between marital states

bysort hhidpn (year): ///
    gen transition = ///
    marital_r != marital_r[_n - 1]

bysort hhidpn (year): ///
    replace transition = 0 if _n == 1


* Identify individuals with no observed transition

bysort hhidpn: ///
    egen NT = ///
    min(cond(transition == 1, 0, 1))


* Count total marital transitions

bysort hhidpn: ///
    egen num_transition = ///
    total(transition)


* Identify transitions from divorce to marriage

bysort hhidpn (year): ///
    gen remarried = ///
    num_transition >= 1 & ///
    transition == 1 & ///
    marital_r[_n - 1] == 2 & ///
    marital_r == 1


* Retain remarriage status during subsequent married waves

bysort hhidpn (year): ///
    replace remarried = 1 ///
    if remarried[_n - 1] == 1 & ///
       remarried == 0 & ///
       marital_r == 1


* Indicator for any observed remarriage experience

bysort hhidpn: ///
    egen remarrExp = ///
    max(remarried == 1)
Code clarification: In this page, marital_s = 1 represents married and marital_s = 2 represents divorced. Therefore, the remarriage transition is coded as a movement from status 2 to status 1.

Section 05

Event-Study Estimation and Visualization

The event-study specification estimates food-insecurity outcomes at several periods before and after divorce. The coefficients are then plotted with 95% confidence intervals to show the timing and persistence of estimated changes.

Event-study regression and coefficient plot Stata 
* Estimate the event-study model

areg foodInsec ///
    tm3 tm2 o.tm1 tm0 tp1 tp2 tp3 ///
    $x yy* ///
    [pweight = baseline_weight] ///
    if divExp == 1, ///
    absorb(hhidpn)

estimates store foodInsec


* Plot event-study coefficients

coefplot ///
    (foodInsec, label("Food Insecurity")), ///
    omitted ///
    drop(_cons $x yy*) ///
    vertical ///
    baselevels ///
    xlabel(, angle(0)) ///
    connect(l) ///
    rename( ///
        tm3 = "-6 or more" ///
        tm2 = "-4" ///
        tm1 = "-2" ///
        tm0 = "0" ///
        tp1 = "+2" ///
        tp2 = "+4" ///
        tp3 = "+6 or more") ///
    ci(95) ///
    xline(4) ///
    yline(0, lpattern(dash)) ///
    graphregion( ///
        fcolor(white) ///
        ifcolor(white) ///
        ilcolor(white)) ///
    ylabel( ///
        -0.1(0.1)0.2, ///
        format(%3.1f) ///
        noticks ///
        nogrid) ///
    ytitle("Change in prevalence of food insecurity") ///
    xtitle("Years before/after gray divorce") ///
    title("Food Insecurity")
Interpretation: The omitted pre-divorce period serves as the reference category. Pre-treatment coefficients can be used to assess trends before divorce, while post-divorce coefficients illustrate the evolution of food insecurity following the marital transition.

Research Code Notice

These examples are simplified and reformatted for demonstration purposes. Variable definitions, sample restrictions, survey weights, identification assumptions, and diagnostic tests should be reviewed carefully before using the code in another research project. The underlying restricted research dataset is not distributed through this website.