** 	THIS FILE COMPARES OBSERVED WITH PREDICTED LAND USE
** 	marea = model area ("predicted")
**	sarea = survey area ("observed")

**	NEDUMARAN
**	15 OCT. 2008

	clear
	set mem 200m

** 	CHANGE THE DIRECTORY
	cd $DIR

** 	OPEN SOLUTION.DTA
	use BSL_Solution, clear

** 	VALIDATE ONLY THE BASELINE SCENARIO
	keep if run == "BSL"

** 	KEEP ONLY THE FIRST 3 YEARS
	keep if period==2&3

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


** 	GROUP THE CROPS

**	RICE AREA
	egen marea_irice = rowtotal(grow_irice_*)
	egen marea_rfrice = rowtotal(grow_rfrice_*)
	egen marea_sprice = rowtotal(grow_sprice_*)
**	TOTAL RICE AREA

	egen marea_rice = rowtotal(grow_irice_* grow_rfrice_* grow_sprice_*)	


**	RAINFED CROPS
	egen marea_rfmaize = rowtotal(grow_rfmaize)
	egen area_mil = rowtotal( grow_cp_mil_1 grow_cp_mil_2 grow_cp_mil_gc_1 grow_cp_mil_gc_2 grow_cp_gc_mil_1 grow_cp_gc_mil_2 grow_millet_wocash) /*Millet and guinea corn*/
	egen area_mil_bean = rowtotal(grow_b_mb_* grow_cp_mil_b_* grow_cp_gmil_b_* grow_cp_gmil2_b_* grow_b_mb_wocash)/*millet and beans*/
	egen area_gnut_bean = rowtotal(grow_b_gb grow_b_gb_m grow_b_gb_fm grow_cp_gnut_b grow_cp_gnut_b_m grow_cp_gnut_b_fm grow_cp_gnut_bb)/*groundnut and beans*/
	egen area_gnut = rowtotal(grow_b_gnut_* grow_cp_gnut) /*Groundnut sole crop*/

**	IRRIGATED CROP AREA
	egen marea_tomato = rowtotal(grow_tom_*)
	egen marea_onion = rowtotal(grow_onion_*)
	egen marea_irmaize = rowtotal(grow_feb_irrimaize_* grow_mar_irrimaize_* grow_aprl_irrimaize_*)

**	ANALYZE THE INTERCROPS

	g ic_mb = area_mil_bean*0.9
	g ic_bm = area_mil_bean*0.1

	g ic_gb = area_gnut_bean*0.9
	g ic_bg = area_gnut_bean*0.1

	g marea_millet = area_mil+ ic_mb
	g marea_beans = ic_bm+ic_bg
	g marea_groundnut = area_gnut+ic_gb


	keep marea_*  clustID period

** 	CALCULATE THE AVERAGE OVER THE PERIODS
	collapse (mean) marea_* ,by (clustID) fast

** 	RESHAPE TO CROPS
	reshape long marea_ ,i(clustID) j(crop) string
	rensfix _

** 	COLLAPSE TO THE CROP LEVEL
	collapse (sum) marea, by(crop) fast
	
**	RENAME CROPS
	g str18 cropname = "."
	replace cropname = "Rice" 			if crop == "rice"
	replace cropname = "Irrigated Rice" 	if crop == "irice"
	replace cropname = "Rainfed Rice" 		if crop == "rfrice"
	replace cropname = "Beans" 			if crop == "beans"
	replace cropname = "Groundnut" 		if crop == "groundnut"
	replace cropname = "Irrigated Maize" 		if crop == "irmaize"
	replace cropname = "Onion" 		if crop == "onion"
	replace cropname = "Tomato" 		if crop == "tomato"
	replace cropname = "Sup Irrigated Rice" 		if crop == "sprice"
	replace cropname = "Rainfed Maize" 		if crop == "rfmaize"
	
	replace cropname = "Millet" 		if crop == "millet"
	replace cropname = "Beans" 		if crop == "beans"


** 	MERGE WITH THE AREAS FROM THE SURVEY
	sort  crop
	merge crop using survey_croparea

	drop _m 

   drop if sarea==.
	drop if marea==.
**	CALCULATE PERCENTAGES 

	local type  "s m"
	
		foreach t of local type {
			quietly egen t_`t'area = total(`t'area) if marea != .
			quietly gen p_`t'area = `t'area / t_`t'area * 100 if marea != .
			drop t_`t'area 
			}


	regress p_marea p_sarea, r noc
	g b  = _b[p_sarea]    
	g se = _se[p_sarea]
	g t = b / se
	g df = e(df_r)
	g sig = tprob(df, t)


** SCATTER: OBSERVED VS. SIMULATED LAND USE (%)
	graph twoway scatter  p_marea p_sarea, msize(medium) mlabel(cropname) ms(square) mlc(black) mc(none) || function y = _b[p_sarea]*x, range(0 40)   ylabel(0(20)60, angle(horizontal) labsize(medium) grid glc(black) glw(*0.5) gmax) xlabel(0(20)80, labsize(medium) grid glc(black) glw(*0.5) gmax) xtitle("Observed values (%)", size(medlarge) margin(medium)) ytitle("Simulated values (%)", size(medlarge) margin(medium)) legend(off) plotregion(fc(white) lc(black)) graphregion(fc(white) margin(b-3 l-5)) ysize(20) xsize(20) title("Observed vs. simulated land use" "as a % of the total area under crops", size(medlarge)) 
	graph export "D:\MPMAS\Ghana_casestudy\TONO_10\Stata\results\Valarea_per_BSL.tif", as(tif) replace



** BAR: OBSERVED AND SIMULATED LAND USE (HA)
	**graph hbar (asis) sarea marea , over(cropname, gap(*0.6) sort(sarea_ws) ) bar(1, fcolor(gs7) lcolor(gs7)) bar(2, fcolor(black) lcolor(black))  plotregion(lc(black)) graphregion(fc(white) margin(b-5 l+0)) ytitle("Crop area (ha)", margin(small)) bargap(0) outergap(*0) ylabel(0(50)350 , glc(black*0.5) glw(*0.5)) ytick(0(25)350) ysize(20) xsize(17) legend(label(1 "Observed") label(2 "Simulated") ring(0) pos(2) col(1) size(medsmall)  nobox symx(*0.5) region(lc(white)))
	**graph export "D:\MPMAS\Ghana_casestudy\TONO_01\stata/results/ValArea_barha_BSL.tif", as(tif) replace

** BAR: OBSERVED AND SIMULATED LAND USE (%)
	**graph hbar (asis) p_marea p_sarea , over(cropname, gap(*0.6) sort(sarea_ws) ) bar(1, fcolor(gs7) lcolor(gs7))  bar(2, fcolor(black) lcolor(black)) plotregion(lc(black)) graphregion(fc(white) margin(b-5 l+0)) yscale(range(0 30)) ytitle("Crop area (%)", margin(small)) bargap(0) outergap(*0.5) ylabel(0(5)60, glc(black*0.5) glw(*0.5) gmax gmin) ymtick(0(1)60) ysize(20) xsize(17) legend(label(1 "Simulated") label(2 "Observed") ring(0) pos(2) col(1) size(medsmall) width() nobox symx(*0.5) region(lc(white)))
	**graph export "D:\MPMAS\Ghana_casestudy\TONO_01\stata/results/ValArea_barp_BSL.tif", as(tif) replace

** BAR: OBSERVED AND SIMULATED LAND USE (%)
	graph bar  p_marea p_sarea , over(cropname, gap(*0.6) sort(sarea_ws) label(angle(forty_five)) ) bar(1, fcolor(forest_green*0.7) lcolor(forest_green*0.7))  bar(2, fcolor(black) lcolor(black)) plotregion(lc(black)) graphregion(fc(white) margin(b-2 l+0)) yscale(range(0 30)) ytitle("Crop area (%)", margin(small)) bargap(0) outergap(*1) ylabel(0(10)60, glc(black*0.5) glw(*0.5) gmax gmin)  ysize(10) xsize(12) legend(label(1 "Simulated") label(2 "Observed") ring(0) pos(2) col(1) size(medsmall) width() nobox symx(*0.5) region(lc(white)))
	graph export "D:\MPMAS\Ghana_casestudy\TONO_10\stata/results/ValArea_barp_BSL.tif", as(tif) replace
gg
