****************************************************************************************
* THIS DO FILE RUNS 2SLS REGRESSIONS TO TEST THE FOOD BUDGET SQUEEZE 
* HYPOTHESIS. THE MAIN REGRESSION USES A POOLED CROSS SECTION DATA SET
* BUILT FROM 4 ROUNDS OF THE NSSO CES (1987, 1993, 2004, 2009).
* DEP VAR: Log CALORIE INTAKE
* KEY INDEP VAR: LOG REAL EXPENDITURE ON FUEL
* INSTRUMENT FOR KEY INDEP VAR: DUMMY VARIABLE FOR MARKET SOURCE OF COOKING ENERGY
* CONTROLS: DIVERSIFICATION, PRICE RATIOS, ACCESS TO SAFE WATER, HOME GROWN CEREALS
*           MEALS EATEN OUTSIDE HOME
* ADDITIONAL CONTROLS: AGE, AGE-SQ, EDUCATION, CASTE, RELIGION, ADJUSTED HOUSEHOLD SIZE
****************************************************************************************


version 13
clear all
set more 1


* ----- USE DATA
use "CES43_50_61_66_pooled_RURAL.dta"

* Labels
label variable lcal_pc "Calorie Intake (kcal per capita)"
label variable cal_pcpd "Calorie Intake (kcal per capita per day)"
label variable rtepc "Total Expenditure (1983 rupees per capita)"
label variable food_all_ratio "Price Ratio (Food/All)"
label variable fdi "Diet Diversification Index"
label variable cer_ncer_ratio_hh "Price Ratio (Cereals/Noncereals)"
label variable safe_water "Access to Safe Water"
label variable hmgrn_cereals "Cereals grown at home (proportion of total)"
label variable hhsize_adj "Household size (adjusted for age and sex)"
label variable head_age "Age of household head (years)"
label variable meals_out "Meals Eaten Outside Home"
label variable cooking_dum "Cooking Source Dummy (Mkt=1, NonMkt=0)"
label variable lrfuel_market_pc "Log Real Fuel Expenditure"
label variable lcer_ncer_ratio_hh "Log Price Ratio (Cereals/Noncereals)" 
label variable lfood_all_ratio "Log Price Ratio (Food/All)"


*----- Drop Outliers 

* (1) in terms of calorie intake
* (Note: drop 1% at the top and bottom of the distribution for all years combined)
summ cal_pcpd [fweight=wt_int], detail
keep if (cal_pcpd>r(p1) & cal_pcpd<r(p99)) 


****************************************************************************************
* ----------------------------- FIGURES ------------------------------------------------

preserve

collapse fe_share nfe_share rfepc rnfepc med_nfe_share edu_nfe_share  ///
  conveyance_nfe_share fuel_nfe_share rmed_pc redu_pc rconveyance_pc rfuel_market_pc ///
   [fweight=wt_int], by(year)


* ------------------------------------ FIGURE 1 ------------------------------- *

* Figure 1(a): Real expenditure (Nonfood + Food)
tw (connected  rfepc year, msymbol(S)) ///
   (connected  rnfepc year, msymbol(T)), graphregion(color(white)) ///
	   legend(order(1 "Food" 2 "Nonfood")) ///
   ytitle("1960-61 rupees") xtitle("")
graph export "foodreal.pdf", replace

* Figure 1(b): Real expenditure (Nonfood)
tw (connected  rmed_pc year, msymbol(O)) ///
   (connected  redu_pc year, msymbol(S)) ///
   (connected  rconveyance_pc year, msymbol(D)) ///
   (connected  rfuel_market_pc year, msymbol(T)), graphregion(color(white)) ///
   legend(order(1 "Health" 2 "Education" 3 "Conveyance" 4 "Fuel")) ///
  ytitle("1960-61 rupees") xtitle("") 
graph export "essentialreal.pdf", replace

	   
* ----------------------------- FIGURE 2 -------------------------------------- *

* Figure 2(a): Share of Food and Nonfood Expenditure 
twoway (connected fe_share year, msymbol(D)) ///
       (connected nfe_share year, msymbol(T) yaxis(2)), graphregion(color(white)) ///
       xtitle("") ///
	   legend(order(1 "Food" 2 "Non-food")) ///
       ytitle(food expenditure (% of MPCE)) ///
       ytitle(nonfood expenditure (% of MPCE), axis(2)) 
graph export "foodshare.pdf", replace

* Figure 2(b): Shares of Nonfood Expenditure
tw (connected  med_nfe_share year, msymbol(O)) ///
   (connected  edu_nfe_share year, msymbol(S)) ///
   (connected  conveyance_nfe_share year, msymbol(D)) ///
   (connected  fuel_nfe_share year, msymbol(T)), graphregion(color(white)) ///
   legend(order(1 "Health" 2 "Education" 3 "Conveyance" 4 "Fuel")) ///
   ytitle("proportion of total expenditure (%)") xtitle("") 
graph export "essentialsharenf.pdf", replace


restore



* --------------------------- FIGURE 3 ---------------------------------------- *

graph bar (mean) cooking_dum [fweight=wt_int], over(year) graphregion(color(white)) ///
 ytitle("mkt source (share of hhlds)") title("Over Years")
graph save cook_year, replace

graph bar (mean) cooking_dum [fweight=wt_int], over(lrtepc_deciles) graphregion(color(white)) ///
 ytitle("mkt source (share of hhlds)") title("Across Expenditure Deciles")
graph save cook_dec, replace

graph combine "cook_year.gph" "cook_dec.gph", graphregion(color(white)) col(2) ///
 title("Use of Market Sources of Cooking Energy")
graph export "energy_source_1.pdf", replace



* --------------------------- FIGURE 4 ---------------------------------------- *
graph bar (mean) rfuel_market_pc [fweight=wt_int] if year==1987, over(cooking_dum) graphregion(color(white)) ///
 ytitle("1960-61 rupees") title("1987-88")
graph save "enrg1.gph", replace

graph bar (mean) rfuel_market_pc [fweight=wt_int] if year==1993, over(cooking_dum) graphregion(color(white)) ///
 ytitle("1960-61 rupees") title("1993-94")
graph save "enrg2.gph", replace

graph bar (mean) rfuel_market_pc [fweight=wt_int] if year==2004, over(cooking_dum) graphregion(color(white)) ///
 ytitle("1960-61 rupees") title("2004-05")
graph save "enrg3.gph", replace

graph bar (mean) rfuel_market_pc [fweight=wt_int] if year==2009, over(cooking_dum) graphregion(color(white)) ///
 ytitle("1960-61 rupees") title("2009-10")
graph save "enrg4.gph", replace

graph combine "enrg1.gph" "enrg2.gph" "enrg3.gph" "enrg4.gph", graphregion(color(white)) col(2) ///
  title("Real Per Capita Expenditure on Fuel") subtitle("(0=Non-Market, 1=Market)")
graph export "energy_source_2.pdf", replace





****************************************************************************************
* -------------------------------------- TABLES ----------------------------------------



* ------------------------------ TABLE 1--------------------------------------- *

eststo clear
qui eststo d1: estpost tabstat cal_pcpd rtepc rfuel_market_pc hhsize_adj head_age meals_out fdi safe_water cooking_dum food_all_ratio  cer_ncer_ratio_hh [fweight=wt_int] if year==1987, listwise statistics(mean sd) columns(statistics)
qui eststo d2: estpost tabstat cal_pcpd rtepc rfuel_market_pc hhsize_adj head_age meals_out fdi safe_water cooking_dum food_all_ratio  cer_ncer_ratio_hh [fweight=wt_int] if year==1993, listwise statistics(mean sd) columns(statistics)
qui eststo d3: estpost tabstat cal_pcpd rtepc rfuel_market_pc hhsize_adj head_age meals_out fdi safe_water cooking_dum food_all_ratio  cer_ncer_ratio_hh [fweight=wt_int] if year==2004, listwise statistics(mean sd) columns(statistics)
qui eststo d4: estpost tabstat cal_pcpd rtepc rfuel_market_pc hhsize_adj head_age meals_out fdi safe_water cooking_dum food_all_ratio  cer_ncer_ratio_hh [fweight=wt_int] if year==2009, listwise statistics(mean sd) columns(statistics)
 
esttab d1 d2 d3 d4, ///
 cells("mean(fmt(3))" "sd(fmt(3) par)") nomtitle nonumber label 


* --------------------------------- TABLE 2  ---------------------------------- * 
* ---- Use the following codes for "cooking": 
* coke, coal-01, firewood and chips-02, LPG-03, gobar gas-04, dung cake-05, charcoal-06, 
* kerosene-07, electricity, others and no cooking arrangement-08

eststo clear
qui eststo a1: estpost tabulate cooking [fweight=wt_int] if year==1987  
qui eststo a2: estpost tabulate cooking [fweight=wt_int] if year==1993  
qui eststo a3: estpost tabulate cooking [fweight=wt_int] if year==2004  
qui eststo a4: estpost tabulate cooking [fweight=wt_int] if year==2009  

* ------------------- other than last row --------------------------- * 
esttab a1 a2 a3 a4 , cells("pct(fmt(2))") ///
 label nonumber nomtitle noobs
 
* ------------------------- last row -------------------------------- * 
* -------- Test for the Equality of Mean Fuel Expenditure by Energy Source (commercial versus noncommercial)
eststo clear
qui eststo e1: estpost ttest rfuel_market_pc if year==1987, by(cooking_dum) unequal
qui eststo e2: estpost ttest rfuel_market_pc if year==1993, by(cooking_dum) unequal
qui eststo e3: estpost ttest rfuel_market_pc if year==2004, by(cooking_dum) unequal
qui eststo e4: estpost ttest rfuel_market_pc if year==2009, by(cooking_dum) unequal
* Manually copy results
esttab e1 e2 e3 e4



* ----------------------- REGRESSION RESULTS 


global controls_new "hhsize_adj head_age head_age2 meals_out caste_dum* religion_dum* edu_dum*"



* -------------------------------------------------------------------------------------------- *
* --------- Within Expenditure Deciles (Within: expenditure deciles, state-region, occupation)

* Drop households that use firewood & chips as primary source of cooking fuel
* Note: this is to address the concern that the cooking dummy might impact calorie 
*       intake through calorie needs

* 10th expenditure decile is the excluded expenditure decile dummy 
drop lrtepc_dec_dum10


* ----------------------------- TABLE 3 ------------------------------------- * 

* ---- Reduced Form and First Stage
* Reduced Form
qui eststo b1: reg lcal_pc cooking_dum lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water year hmgrn_cereals $controls_new [fweight=wt_int], vce(cluster srcode3866)

* First Stage Regression: Log real fuel_market
qui eststo b2: reg lrfuel_market_pc cooking_dum lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water year hmgrn_cereals $controls_new [fweight=wt_int], vce(cluster srcode3866)

esttab b1 b2, keep(cooking_dum lcer_ncer_ratio_hh lfood_all_ratio fdi hhsize_adj head_age head_age2 year) ///
cells(b(star fmt(%9.3f)) se(par fmt(%9.3f)) ) ///
starlevels(* 0.10 ** 0.05 *** 0.01) style(tex) label ///
stats(r2_a N, fmt(2 0) label("Adjusted R-sq" "Observations")) nonotes


* ----------------------------- TABLE 4 ------------------------------------- * 

* Basic model with diversification index
qui eststo a1: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* $controls_new [fweight=wt_int], cluster(srcode3866) partial(st_rgn*)

* Add occupation dummies
qui eststo a2: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* $controls_new [fweight=wt_int], cluster(srcode3866) partial(st_rgn* occp*)

* Add safe water
qui eststo a3: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water $controls_new [fweight=wt_int], cluster(srcode3866) partial(st_rgn* occp*)

* Add home grown cereals
qui eststo a4: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water hmgrn_cereals $controls_new [fweight=wt_int], cluster(srcode3866) partial(st_rgn* occp*)

* Add time trend
qui eststo a5: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water hmgrn_cereals $controls_new year [fweight=wt_int], cluster(srcode3866) partial(st_rgn* occp*)

esttab a1 a2 a3 a4 a5 , keep(lrfuel_market_pc lcer_ncer_ratio_hh lfood_all_ratio  fdi hhsize_adj head_age head_age2) ///
cells(b(star fmt(%9.3f)) se(par fmt(%9.3f)) ) label ///
starlevels(* 0.10 ** 0.05 *** 0.01) style(tex) nomtitle ///
stats(r2_a idp N, fmt(3 3 3 0) label("Adjusted R-sq" "Underid p-val" "Observations")) nonotes


* ----------------------------- TABLE 5 ------------------------------------- * 

* --- Robustness checks
preserve
drop if cooking==2
qui eststo r1: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water hmgrn_cereals $controls_new year [fweight=wt_int], cluster(srcode3866) partial(st_rgn* occp*)
restore

preserve
drop if (fuelwood_free==1|fuelwood_free==.)
qui eststo r2: ivreg2 lcal_pc (lrfuel_market_pc = cooking_dum) lrtepc_dec_dum* fdi lcer_ncer_ratio_hh lfood_all_ratio st_rgn* occp* safe_water hmgrn_cereals $controls_new year [fweight=wt_int], cluster(srcode3866) partial(st_rgn* occp*)
restore

esttab r1 r2 , keep(lrfuel_market_pc lcer_ncer_ratio_hh lfood_all_ratio  fdi hhsize_adj head_age head_age2) ///
cells(b(star fmt(%9.3f)) se(par fmt(%9.3f)) ) label ///
starlevels(* 0.10 ** 0.05 *** 0.01) style(tex) nomtitle ///
stats(r2_a idp N, fmt(3 3 3 0) label("Adjusted R-sq" "Underid p-val" "Observations")) nonotes