/** Replicate mortality ratio calculation of Weir (2016) for a given year
This only calculates the numerator and denominator summed by age-demographic group
*/


global year : env YEAR

* specify which death date variables to use: rai (HRS), ran (NDI), or rad (NDI or HRS if NDI missing)
local dthdtvar rad

log using logs/weir_mortality_repl_yr$year.log, replace

* Compute exposure time at each age for this year in HRS
use data/hrs_analytic_P.dta, clear

** drop the already-dead to speed up processing
drop if !missing(`dthdtvar'year) & `dthdtvar'year < $year

** drop the not-yet-eligible to speed up processing
*** Weir: "begin observation on each individual with her first interview taken after her birth cohort became age-eligible"
drop if missing(insampyr) | insampyr > $year

** compute age at start of year
gen age = min($year - rabyear - 1, 100) if !missing(rabyear)

** expand each person to get an exposure observation for each age pre- and post- birthday during the year
expand 2, gen(postbday)
replace age = min(age + postbday, 100)

** drop observations at age 100 and older
drop if age >= 100

** flag ages when people actually died
gen died = `dthdtvar'year==$year & !missing(`dthdtvar'month) & ((postbday==0 & `dthdtvar'month <= rabmonth) | (postbday==1 & `dthdtvar'month >= rabmonth))
*** split deaths between ages if death occurs in same month as birthday
replace died = died/2.0 if `dthdtvar'year==$year & !missing(`dthdtvar'month) & `dthdtvar'month==rabmonth

** compute exposure at each age during the year
gen exp_start = 1 if postbday==0
gen exp_end = rabmonth if postbday==0

replace exp_start = rabmonth if postbday==1
replace exp_end = 12 if postbday==1

replace exp_start = max(exp_start, insampmo) if insampyr==$year
replace exp_end = min(exp_end, `dthdtvar'month) if `dthdtvar'year==$year & !missing(`dthdtvar'month)

gen exposure = max(exp_end - exp_start + 0.5, 0)/12.0

* Use life table to compute expected number of deaths at individual level
sort black male age
preserve
use ../../base_data/cdc_nvsr_lifetable$year.dta, clear
cap confirm var hispanic
if _rc == 0 {
   drop if !missing(hispanic)
}
tempfile lifetab$year
save `lifetab$year'
restore

merge m:1 black male age using `lifetab$year'
drop if _merge==2

li if _merge==1

** Q: should people with missing gender/race be dropped?
* drop if _merge==1

* generate morality rate m_x
gen mx = dx / Lx

gen expdied = exposure * mx

li if _n <= 250

* Aggregate actual and expected number of deaths across individuals by demographic group
collapse (sum) exposure died expdied dx Lx, by(black male age)
gen year=$year
li

save data/weir_mortality_repl$year.dta, replace