*******************************************************************************
*To accompany Knittel and Metaxoglou (2008)
*Work in Progess. Please do not cite without Permisssion
*Knittel      : crknittel@ucdavis.edu
*Metaxoglou   : konstantinos.metaxoglou@bateswhite.com
*******************************************************************************
clear
set memo 10000m
set more off
capture log close

*******************************************************************************
*Define globals for directories and files
*******************************************************************************
global path_main          	"H:\RCOptim\Results\~BLP delivery CK\~BLP delivery"
global path_temp          	"e:\crknittel"
global path_in            	"$path_main"
global path_out	  	  	"$path_main"

global file_in_optim            "$path_in\Optimization results\Optimization results.dta"
global file_optim_tmp		  "$path_temp\blp_optim_results_tmp.dta"

global file_in_mkt_power        "$path_in\market power results\blp_mkt_power_results_all.dta"
global file_in_mkt_power_tmp    "$path_temp\blp_mkt_power_results_all_tmp.dta"

global file_in_merger           "$path_in\merger results\blp_merger_results_all.dta"
global file_in_merger_tmp       "$path_temp\blp_merger_results_all_tmp.dta"

global file_BLP_analysis        "$path_temp\BLP_analysis.dta"
global file_BLP_analysis_final  "$path_temp\BLP_analysis_final.dta"

global path_tables              "H:\RCOptim\Results\BLP_base\Tables"
global path_figures              "$path_main\figures\"
global figure17                  "$path_temp\figure17.dta"
global local_minima              "$path_temp\blp_local minima.dta"

*******************************************************************************
*Load market power results - BLP
*******************************************************************************
use "$file_in_mkt_power", clear
capture rename brand product
sort optmethod stvalue market product
save "$file_in_mkt_power_tmp", replace

*******************************************************************************
*Load merger results - BLP
*******************************************************************************
use "$file_in_merger", clear
keep if optmethod~=6

capture rename start stvalue
capture rename model_id product
capture drop *chng

sort optmethod stvalue market product
save "$file_in_merger_tmp", replace

*******************************************************************************
*Merge the optimization, convergence, market power and merger results - BLP
*Prepare a dataset to produce figures and tables
*******************************************************************************

use "$file_in_optim", clear

*******************************************************************************
*Convergence 
*******************************************************************************
capture gen conv2=0
replace conv2=1 if optmethod==1     & exitinfo==1
replace conv2=1 if optmethod==1     & exitinfo==2
replace conv2=1 if optmethod==2     & exitinfo==1
replace conv2=1 if optmethod==3     & exitinfo<4000
replace conv2=1 if optmethod==4     & exitinfo==2
replace conv2=1 if optmethod==5     & exitinfo==2
replace conv2=1 if optmethod==5     & exitinfo==4
replace conv2=1 if optmethod==7  	& fcnevals<4000
replace conv2=1 if optmethod==8     & exitinfo==1
replace conv2=1 if optmethod==9     & exitinfo==1
replace conv2=1 if optmethod==10    & exitinfo==1	& fcnevals==1020

replace conv2=1 if stvalue>2		& stvalue<17	& optmethod==2
replace conv2=1 if stvalue==21	& optmethod==2
replace conv2=1 if stvalue==22	& optmethod==2
replace conv2=1 if stvalue==24	& optmethod==2
replace conv2=1 if stvalue==25	& optmethod==2
replace conv2=1 if stvalue==26	& optmethod==2
replace conv2=1 if stvalue==28	& optmethod==2
replace conv2=1 if stvalue==30	& optmethod==2
replace conv2=1 if stvalue==31	& optmethod==2
replace conv2=1 if stvalue==33	& optmethod==2
replace conv2=1 if stvalue==34	& optmethod==2
replace conv2=1 if stvalue==36	& optmethod==2
replace conv2=1 if stvalue==37	& optmethod==2
replace conv2=1 if stvalue>38		& stvalue<48	& optmethod==2
replace conv2=1 if stvalue==49	& optmethod==2

save "$file_optim_tmp", replace
*******************************************************************************

sort  optmethod stvalue

merge optmethod stvalue using "$file_in_mkt_power_tmp"
tab   _merge
keep if _merge==3
drop _merge

sort  optmethod stvalue market product 
merge optmethod stvalue market product using "$file_in_merger_tmp"
tab _merge
keep if _merge==3
drop _merge
compress
save "$file_BLP_analysis", replace

use "$file_BLP_analysis", clear

gen     flag=0

*******************************************************************************
*Identify products with market shares
*corresponding to the four quartiles of the market share distribution
*******************************************************************************
sort optmethod stvalue share_obs market product, stable

global  market_25=market  in 554
global product_25=product in 554

global  market_50=market  in 1109
global product_50=product in 1109

global  market_75=market  in 1663
global product_75=product in 1663

global  market_100=market  in 2217
global product_100=product in 2217

gen flag_025  = (market==$market_25) *(product==$product_25)
gen flag_050  = (market==$market_50) *(product==$product_50)
gen flag_075  = (market==$market_75) *(product==$product_75)
gen flag_100  = (market==$market_100)*(product==$product_100)

table market product if flag_025 ==1, c(min share_obs) f(%8.4fc)
table market product if flag_050 ==1, c(min share_obs) f(%8.4fc)
table market product if flag_075 ==1, c(min share_obs) f(%8.4fc)
table market product if flag_100 ==1, c(min share_obs) f(%8.4fc)

*******************************************************************************
*Identify the results that correspond to the minimum value
*of the objective function for each of the 10 optimization routines
*******************************************************************************
gsort optmethod fval stvalue -conv2

egen fval_rank=group(optmethod fval stvalue conv2)
egen fval_rank_min=min(fval_rank), by(optmethod conv2)
replace fval_rank=fval_rank-fval_rank_min+1

gen      fval_optmethod_min_flag = 0
replace  fval_optmethod_min_flag = 1 if fval_rank==1&conv2==1

sum      fval if fval_optmethod_min_flag

*******************************************************************************
*Identify the results that correspond to the minimum value
*of the objective function for each of the 10 optimization routines
*for the first 20 starting values
*******************************************************************************

gen stvalue_over_20=stvalue>20
 
gsort optmethod stvalue_over_20 fval stvalue

egen fval_rank_20=group(optmethod fval stvalue stvalue_over_20)
egen fval_rank_min_20=min(fval_rank_20), by(optmethod stvalue_over_20)
gen fval_optmethod_min_flag_20=0
replace fval_optmethod_min_flag_20=1 if fval_rank_20==fval_rank_min_20&stvalue_over_20==0

sum      fval stvalue if fval_optmethod_min_flag_20==1

*******************************************************************************
*Identify the closest substitute for the products with market shares
*corresponding to the four quartiles of the market share distribution
*******************************************************************************
aorder
global elast_cross_all elast_cross1-elast_cross150

*the presence of missing values affects the max calculation
forvalues i=1(1)150 {
	qui replace elast_cross`i'=0 if elast_cross`i'==.	
}

qui sum fval, detail, if conv2==1
qui return list
global fval_upp=r(p90)

*sum elast_cross1-elast_cross150 if flag_025==1  &conv2==1  &fval<$fval_upp
*sum elast_cross1-elast_cross150 if flag_050==1  &conv2==1  &fval<$fval_upp
*sum elast_cross1-elast_cross150 if flag_075==1  &conv2==1  &fval<$fval_upp
*sum elast_cross1-elast_cross150 if flag_100==1  &conv2==1  &fval<$fval_upp
***********HARD CODED************
global cross025=42
global cross050=31
global cross075=125
global cross100=21

gen price_chng  = price_post-price_pre

compress
sort optmethod stvalue market product
save "$file_BLP_analysis_final",replace

*******************************************************************************
use "$file_BLP_analysis_final",clear

keep mean_CV conv2 optmethod stvalue fval market fval_optmethod_min_flag
keep if market>10
collapse (mean) mean_CV fval_optmethod_min_flag (max) conv2 fval, by(optmethod stvalue)

gen 	mean_CV_pop=mean_CV*260

qui sum mean_CV_pop, detail, if conv2==1
qui return list
global low=r(p1)
global upp=0

gen flag=(conv2==1)*(mean_CV_pop>$low)*(mean_CV_pop<$upp)

save "$figure17",replace

*******************************************************************************

use "$file_BLP_analysis_final",clear
*egen gradients1_norm_max=rowmax(gradients1_*)
*egen gradients1_norm_min=rowmin(gradients1_*)
*replace gradients1_norm_min=abs(gradients1_norm_min)
*egen gradients1_norm_inf =rowmax(gradients1_norm_min gradients1_norm_max)

keep gradients1_norm_inf hessians_* conv2 optmethod stvalue fval

collapse (max) gradients1_norm_inf hessians* conv2 fval, by(optmethod stvalue)
egen hessians_min=rowmin(hessians_*)
gen flag_30=(gradients1_norm_inf<=30)*(hessians_min>=0)*conv2
keep if flag_30==1
gsort optmethod stvalue
keep optmethod stvalue flag_30 conv2
compress
save "$local_minima",replace

cd "$path_figures"

use "$file_BLP_analysis_final",clear

*******************************************************************************
*Figure 1 - BLP
*******************************************************************************
qui sum fval, detail, if conv2==1
qui return list
global fval_upp=r(p90)

graph hbox fval, over(optmethod) ytitle(GMM Objective Function Value), if fval<$fval_upp&product==1&market==16
graph export Figure01.wmf ,replace

*******************************************************************************
*Figure 2 - BLP
*******************************************************************************
qui sum fval, detail, if conv2==1
qui return list
global fval_upp=r(p90)
#delimit ;
histogram fval if conv2==1 &fval<$fval_upp&product==1&market==16, 
	freq 
	w(5)
	xtitle(GMM Objective Values); 
graph export Figure02.wmf ,replace;
#delimit cr
*******************************************************************************
*Figure 5 - BLP
*******************************************************************************
global upp=0
qui sum elast_own, detail, if flag_025==1&conv2==1 &fval<$fval_upp
qui return list
global low_025=r(p10)

qui sum elast_own, detail, if flag_050==1&conv2==1 &fval<$fval_upp
qui return list
global low_050=r(p10)

qui sum elast_own, detail, if flag_075==1&conv2==1 &fval<$fval_upp
qui return list
global low_075=r(p10)

qui sum elast_own, detail, if flag_100==1&conv2==1 &fval<$fval_upp
qui return list
global low_100=r(p10)

global w_025=abs($low_025)/25
global w_025=abs($low_050)/25
global w_025=abs($low_075)/25
global w_025=abs($low_100)/25

#delimit ;
histogram elast_own if flag_025==1 & conv2==1 & elast_own>=$low_025  & elast_own<$upp	&fval<$fval_upp, 
	w($w_025)
	title(25th Percentile) 
	xtitle(Own-Price Elasticity) 
	saving(BLP_25, replace) 
	xlabel(,format(%9.2f))
	freq;

histogram elast_own if flag_050==1 & conv2==1 & elast_own>=$low_050  & elast_own<$upp	&fval<$fval_upp, 
	w($w_025)
	title(50th Percentile) 
	xtitle(Own-Price Elasticity) 
	saving(BLP_median, replace) 
	xlabel(,format(%9.2f))
	freq;
	
histogram elast_own if flag_075==1 & conv2==1 & elast_own>=$low_075  & elast_own<$upp	&fval<$fval_upp, 
	w($w_025)
	title(75th Percentile) 
	xtitle(Own-Price Elasticity) 
	saving(BLP_75, replace) 
	xlabel(,format(%9.2f))
	freq;

histogram elast_own if flag_100==1 & conv2==1 & elast_own>=$low_100  & elast_own<$upp	&fval<$fval_upp, 
	w($w_025)
	title(Maximum) 
	xtitle(Own-Price Elasticity) 
	saving(BLP_largest, replace) 
	xlabel(,format(%9.2f))
	freq;

graph combine BLP_25.gph BLP_median.gph BLP_75.gph BLP_largest.gph;
graph export Figure05.wmf,replace;
#delimit cr

erase BLP_25.gph
erase BLP_median.gph 
erase BLP_75.gph
erase BLP_largest.gph

egen own_stddev=sd(elast_own), by(product market conv2)
egen own_stddev_best=sd(elast_own), by(product market conv2 fval_optmethod_min_flag)
sum own_stddev elast_own, detail, if conv2==1
sum own_stddev_best elast_own, detail, if conv2==1&fval_optmethod_min_flag==1

*******************************************************************************
*Figure 6 - BLP
*******************************************************************************
qui sum elast_own, detail, if conv2==1	&fval_optmethod_min_flag==1
qui return list
global low=r(p10)
global upp=0

replace flag=(fval_optmethod_min_flag)*(elast_own>$low)*(elast_own<$upp)*(fval<$fval_upp)

#delimit;
twoway 
(kdensity elast_own                           if optmethod==1  & flag==1)
(kdensity elast_own, lpattern(tight_dot),     if optmethod==2  & flag==1)
(kdensity elast_own, lpattern(dash),          if optmethod==3  & flag==1)
(kdensity elast_own, lpattern(dot),           if optmethod==4  & flag==1)
(kdensity elast_own, lpattern(longdash),      if optmethod==5  & flag==1)
(kdensity elast_own, lpattern(longdash_dot),  if optmethod==6  & flag==1)
(kdensity elast_own, lpattern(dash_3dot),     if optmethod==7  & flag==1)
(kdensity elast_own, lpattern(dash_dot),      if optmethod==8  & flag==1)
(kdensity elast_own, lpattern(shortdash),     if optmethod==9  & flag==1)
(kdensity elast_own, lpattern(shortdash_dot), if optmethod==10 & flag==1),
legend(lab(1 "Quasi-Newton 1") 
       lab(2 "Simplex")
       lab(3 "SolvOpt")
       lab(4 "Conjugate Gradient")
       lab(5 "Quasi-Newton 2")
       lab(6 "JBES GA")
       lab(7 "Simulated Annealing")
       lab(8 "MADS")
       lab(9 "GPS")
       lab(10 "Matlab GA"))  
ytitle(Density) 
xtitle(Own-Price Elasticity);
graph export Figure06.wmf,replace;
#delimit cr

*******************************************************************************
*Figure 7 - BLP
*******************************************************************************
qui sum elast_own, detail, if conv2==1	&fval_optmethod_min_flag_20==1
qui return list
global low=r(p10)
global upp=0

replace flag=(fval_optmethod_min_flag_20)*(elast_own>$low)*(elast_own<$upp)*(fval<$fval_upp)

#delimit;
twoway 
(kdensity elast_own                           if optmethod==1  & flag==1)
(kdensity elast_own, lpattern(tight_dot),     if optmethod==2  & flag==1)
(kdensity elast_own, lpattern(dash),          if optmethod==3  & flag==1)
(kdensity elast_own, lpattern(dot),           if optmethod==4  & flag==1)
(kdensity elast_own, lpattern(longdash),      if optmethod==5  & flag==1)
(kdensity elast_own, lpattern(longdash_dot),  if optmethod==6  & flag==1)
(kdensity elast_own, lpattern(dash_3dot),     if optmethod==7  & flag==1)
(kdensity elast_own, lpattern(dash_dot),      if optmethod==8  & flag==1)
(kdensity elast_own, lpattern(shortdash),     if optmethod==9  & flag==1)
(kdensity elast_own, lpattern(shortdash_dot), if optmethod==10 & flag==1),
legend(lab(1 "Quasi-Newton 1") 
       lab(2 "Simplex")
       lab(3 "SolvOpt")
       lab(4 "Conjugate Gradient")
       lab(5 "Quasi-Newton 2")
       lab(6 "JBES GA")
       lab(7 "Simulated Annealing")
       lab(8 "MADS")
       lab(9 "GPS")
       lab(10 "Matlab GA"))  
ytitle(Density) 
xtitle(Own-Price Elasticity);
graph export Figure07.wmf,replace;
#delimit cr

*******************************************************************************
*Figure 11 - BLP
*******************************************************************************
global low =0
qui sum elast_cross$cross025, detail, if flag_025==1&conv2==1
qui return list
global upp_025=r(p90)

qui sum elast_cross$cross050, detail, if flag_050==1&conv2==1
qui return list
global upp_050=r(p90)

qui sum elast_cross$cross075, detail, if flag_075==1&conv2==1
qui return list
global upp_075=r(p90)

qui sum elast_cross$cross100, detail, if flag_100==1&conv2==1
qui return list
global upp_100=r(p90)
global low_100=r(p10)


#delimit;
histogram elast_cross$cross025 if flag_025==1 & conv2==1 & elast_cross$cross025>=$low & elast_cross$cross025<=$upp_025 &fval<$fval_upp, 
	title(25th Percentile) 
	xtitle(Cross-Price Elasticity for Closest Substitute) 
	saving(BLP_cross_25, replace) 
	xlabel(,format(%9.3f))
	freq;

histogram elast_cross$cross050 if flag_050==1 & conv2==1 & elast_cross$cross050>=$low & elast_cross$cross050<=$upp_050 &fval<$fval_upp, 
	title(50th Percentile) 
	xtitle(Cross-Price Elasticity Closest for Substitute) 
	saving(BLP_cross_median, replace) 
	xlabel(,format(%9.3f))
	freq;

histogram elast_cross$cross075 if flag_075==1 & conv2==1 & elast_cross$cross075>=$low & elast_cross$cross075<=$upp_075 &fval<$fval_upp, 
	title(75th Percentile) 
	xtitle(Cross-Price Elasticity for Closest Substitute) 
	saving(BLP_cross_75, replace) 
	xlabel(,format(%9.3f))
	freq;

histogram elast_cross$cross100 if flag_100==1 & conv2==1 & elast_cross$cross100>=$low_100 & elast_cross$cross100<=$upp_100 &fval<$fval_upp, 
	title(Maximum) 
	xtitle(Cross-Price Elasticity for Closest Substitute) 
	saving(BLP_cross_largest, replace) 
	xlabel(,format(%9.3f))
	freq;

graph combine BLP_cross_25.gph BLP_cross_median.gph BLP_cross_75.gph BLP_cross_largest.gph;

graph export Figure11.wmf,replace;
#delimit cr

erase BLP_cross_25.gph
erase BLP_cross_median.gph 
erase BLP_cross_75.gph
erase BLP_cross_largest.gph

*******************************************************************************
* Figure 17 - BLP
*******************************************************************************
use  "$figure17", clear

histogram mean_CV_pop if flag==1, w(.5) xtitle("Change in Consumer Welfare (billions of dollars)") freq
graph export figure17.wmf,replace 

use "$file_BLP_analysis_final", clear
egen price_chng_stddev=sd(elast_own), by(product market conv2)
egen price_chng_stddev_best=sd(elast_own), by(product market conv2 fval_optmethod_min_flag)
sum price_chng_stddev price_chng, detail, if conv2==1
sum price_chng_stddev_best price_chng, detail, if conv2==1&fval_optmethod_min_flag==1


*******************************************************************************
*Figure 20 - BLP
*******************************************************************************
use "$file_BLP_analysis_final", clear
keep optmethod stvalue elast_own fval product market conv2

sort     optmethod stvalue
merge    optmethod stvalue using "$local_minima"
tab     _merge
keep if _merge==3
drop    _merge

drop if flag_30==0
drop if conv2==0

qui sum elast_own, detail, if conv2==1
qui return list
global low=r(p10)
global upp=0

#delimit;
twoway 
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	2	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	3	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	6	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	15	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	16	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	19	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	21	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	23	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	33	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	37	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	40	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	1	&stvalue==	46	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	2	&stvalue==	6	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	2	&stvalue==	11	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	2	&stvalue==	13	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	2	&stvalue==	25	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	2	&stvalue==	44	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	3	&stvalue==	30	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	4	&stvalue==	5	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	4	&stvalue==	13	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	4	&stvalue==	30	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	4	&stvalue==	31	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	4	&stvalue==	37	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	2	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	3	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	4	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	7	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	10	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	12	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	13	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	16	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	28	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	32	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	33	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	37	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	41	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	43	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	44	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	45	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	47	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	5	&stvalue==	49	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	11	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	14	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	16	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	17	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	19	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	22	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	23	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	24	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	25	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	27	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	28	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	35	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	36	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	37	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	40	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	41	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	44	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	8	&stvalue==	49	& elast_own>$low	&	elast_own<$upp 	)
(kdensity	elast_own	if 	optmethod==	9	&stvalue==	25	& elast_own>$low	&	elast_own<$upp 	)
,  legend(off) ytitle(Density) xtitle(Own-Price Elasticity);
graph export Figure20.wmf,replace;

*******************************************************************************
* Figure 21 - BLP
*******************************************************************************
use  "$figure17", clear

sort   fval optmethod stvalue, stable
global fval_drop=fval in 1

drop if fval==$fval_drop

sort     optmethod stvalue
merge    optmethod stvalue using "$local_minima"
tab     _merge
keep if _merge==3
drop    _merge
drop if conv2==0

drop if flag_30==0

table stvalue optmethod

histogram mean_CV_pop if flag_30==1, w(.1) xtitle("Change in Consumer Welfare (billions of dollars)") freq
graph export Figure21.wmf,replace

*******************************************************************************
*Table 2 BLP parameter estimates for the best of the best
*******************************************************************************
#delimit ;

browse  price_mean  const_mean hpwt_mean air_mean mpd_mean space_mean  price_sigma 
const_sigma hpwt_sigma air_sigma mpd_sigma  price_mean_se  const_mean_se hpwt_mean_se 
air_mean_se mpd_mean_se space_mean_se  price_sigma_se  const_sigma_se hpwt_sigma_se 
air_sigma_se mpd_sigma_se if  fval_optmethod_min_flag==1 &product==16 &market==1;

#delimit cr

*******************************************************************************
*Table 4 BLP 
*******************************************************************************

cd "$path_tables"

foreach i in 025 050 075 100	{ 
	use "$file_BLP_analysis_final",clear	
	qui sum fval, detail, if conv2==1
	qui return list
	global fval_upp=r(p90)
	drop if flag_`i'~=1
	drop if conv2~=1
	drop if fval>$fval_upp
	collapse (min) min_elast=elast_own (max) max_elast=elast_own, by(optmethod)
	outsheet using BLP_`i'_elast.csv, comma replace
	
	use "$file_BLP_analysis_final",clear	
	drop if flag_`i'~=1
	drop if conv2~=1
	drop if fval_optmethod_min_flag~=1	
	drop if fval>$fval_upp
	collapse (min) min_elast=elast_own (max) max_elast=elast_own, by(optmethod)
	outsheet using BLP_`i'_best_elast.csv, comma replace
	}

*******************************************************************************
*Table 6 BLP 
*******************************************************************************

global cross025=42
global cross050=31
global cross075=125
global cross100=21


	use "$file_BLP_analysis_final",clear	
	qui sum fval, detail, if conv2==1
	qui return list
	global fval_upp=r(p90)
	drop if flag_025~=1
	drop if conv2~=1
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross025 (max) max_cross=elast_cross$cross025, by(optmethod)
	outsheet using BLP_025_cross_elast.csv, comma replace
	
	use "$file_BLP_analysis_final",clear	
	drop if flag_025~=1
	drop if conv2~=1
	drop if fval_optmethod_min_flag~=1	
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross025 (max) max_cross=elast_cross$cross025, by(optmethod)
	outsheet using BLP_025_best_cross_elast.csv, comma replace

	use "$file_BLP_analysis_final",clear	
	drop if flag_050~=1
	drop if conv2~=1
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross050 (max) max_cross=elast_cross$cross050, by(optmethod)
	outsheet using BLP_050_cross_elast.csv, comma replace
	
	use "$file_BLP_analysis_final",clear	
	drop if flag_050~=1
	drop if conv2~=1
	drop if fval_optmethod_min_flag~=1	
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross050 (max) max_cross=elast_cross$cross050, by(optmethod)
	outsheet using BLP_050_best_cross_elast.csv, comma replace

	use "$file_BLP_analysis_final",clear	
	drop if flag_075~=1
	drop if conv2~=1
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross075 (max) max_cross=elast_cross$cross075, by(optmethod)
	outsheet using BLP_075_cross_elast.csv, comma replace
	
	use "$file_BLP_analysis_final",clear	
	drop if flag_075~=1
	drop if conv2~=1
	drop if fval_optmethod_min_flag~=1	
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross075 (max) max_cross=elast_cross$cross075, by(optmethod)
	outsheet using BLP_075_best_cross_elast.csv, comma replace

	use "$file_BLP_analysis_final",clear	
	drop if flag_100~=1
	drop if conv2~=1
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross100 (max) max_cross=elast_cross$cross100, by(optmethod)
	outsheet using BLP_100_cross_elast.csv, comma replace
	
	use "$file_BLP_analysis_final",clear	
	drop if flag_100~=1
	drop if conv2~=1
	drop if fval_optmethod_min_flag~=1	
	drop if fval>$fval_upp
	collapse (min) min_cross=elast_cross$cross100 (max) max_cross=elast_cross$cross100, by(optmethod)
	outsheet using BLP_100_best_cross_elast.csv, comma replace

*******************************************************************************
*Generating standard deviation of cross-price estimates
*******************************************************************************
use "$file_BLP_analysis_final", clear
drop if conv2==0
drop if fval>$fval_upp
sete more off
forvalues i=1(1)150	{	
	qui replace elast_cross`i'=. if elast_cross`i'==0
	}
stack elast_cross1 - elast_cross150, into(elast_cross_stack) clear
sum elast_cross_stack, detail, if elast_cross_stack~=.

use "$file_BLP_analysis_final", clear
drop if conv2==0
forvalues i=1(1)150	{	
	qui replace elast_cross`i'=. if elast_cross`i'==0
	egen elast_cross`i'_stdev=sd(elast_cross`i'), by(product market)
	}
stack elast_cross1_stdev - elast_cross150_stdev, into(elast_cross_stdev_stack) clear
sum elast_cross_stdev_stack, detail, if elast_cross_stdev_stack~=.

use "$file_BLP_analysis_final", clear
drop if conv2==0
drop if fval_optmethod_min_flag==0
forvalues i=1(1)150	{	
	qui replace elast_cross`i'=. if elast_cross`i'==0
	}
stack elast_cross1 - elast_cross150, into(elast_cross_best_stack) clear
sum elast_cross_best_stack, detail, if elast_cross_best_stack~=.

use "$file_BLP_analysis_final", clear
drop if fval_optmethod_min_flag==0
forvalues i=1(1)150	{	
	egen elast_cross`i'_stdev_best=sd(elast_cross`i'), by(product market conv2 fval_optmethod_min_flag)
	}
stack elast_cross1_stdev_best - elast_cross150_stdev_best if conv2==1&fval_optmethod_min_flag==1, into(elast_cross_stdev_best_stack)
sum elast_cross_stdev_best_stack, detail