clear 
cap log close
set more off
cd "H:\Charles\Income maps"

*shp2dta using "regions2_83_onworldmapbasis.shp", data(world-d) coor(world-c) genid(id)

	use All_Kenny.dta, clear
		*Make corrections to regions
			replace Region = "TIran, Islamic Rep.a" 	 							if Region=="Tirana"
			replace Region = "Yukon Territory, Northwest Territories, and Nunav*"	if Region=="Yukon Territory, Northwest Territories, and Nunavut"
			replace Country = "Germany" if regexm(Country, "Germany")    
			renvars, lower
			tempfile gdp
			save 	 `gdp'
			
	use world-d.dta, clear
		drop if Country=="Russian Federationn Federation"						//Overwrites other regions in Russia
			replace Region  = regexr(Region, "Bosnia and Herzegovinaar", "Bihar")
			replace Region  = regexr(Region, "Bosnia and Herzegovinaor", "Bihor")
			replace Region  = "Chittagong"	if Region=="Chittagong H.T.Chittagong H.T."
			replace Region  = "Mymensingh"	if Region=="MymensinghMymensingh"
			replace Region  = "Rangpur" 	if Region=="Ranpur"
			replace Region  = "Tangail" 	if Region=="Tangali"
			replace Country = "Germany" if regexm(Country, "Germany")
			renvars, lower
			tempfile world
			save `world'
			
		mmerge country region using `gdp'
		drop if _merge==2 
		drop _merge
		ren regionalgdppc gdp
			
		egen year_last = max(year), by(region)
		drop if year!=year_last
		drop if year_last<2000			//There are 22 observations from 1990. This is a stretch
		
		
		_pctile gdp, p(10 (10) 90)
		return list						//List categories so that they can be repeated with WDI data
		xtile gdpcat = gdp, nq(10)
		tab gdpcat

	 spmap gdpcat using "world-c.dta", id(id)	///
		  fcolor("220 0 20" "208 0 47" "185 0 70"  "162 0 93" "140 0 100" "120 0 120" "112 0 143" "86 0 169" "64 0 191" "47 0 208" )  clm(unique) ///
		  legend(label(2 "< $2000") label(3 "$2000-3500") label(4 "$3500-4700") label(5 "$4700-5900") label(6 "$5900-7700") label(7 "$7700-10,800")  ///
		  label(8 "$10,800-17,600")  label(9 "$17,600-26,200") label(10 "$26,200-33,900") label(11 "> $33,900") ///
		  label(12 "No data") size(vsmall) title("Regional GDP per capita", size(vsmall))) ndfcolor("192 192 192") ndocolor("192 192 192") ///
		  ocolor(none ...) saving(g1, replace) 
			
		renvars  gdp gdpcat /  gdpr gdpcatr
		tempfile region
		save `region'

*WDI data

	import excel "H:\Charles\Income maps\country gdp.xls", sheet("Sheet1") firstrow clear
		renvars, lower
		keep countryname ba
		renvars countryname ba / country gdp
		lab var gdp "GDP per capita 2010"
		
			replace country = "Syria"							if country=="Syrian Arab Republic"
			replace country = "Venezuela"							if country=="Venezuela, RB"
			replace country = "Korea, Rep."							if country=="Republic of Korea"
			replace country = "North Korea"							if country=="Korea, Dem. Rep."
			replace country = "Macedonia"							if country=="Macedonia, FYR"
			replace country = "Micronesia"							if country=="Micronesia, Fed. Sts."
			replace country = "Macao"							if country=="Macao SAR, China"
			replace country = "Laos"							if country=="Lao PDR"
			replace country = "Gambia"							if country=="Gambia, The"
			replace country = "Hong Kong"							if country=="Hong Kong SAR, China"
			replace country = "Democratic Republic of the Congo" 				if country=="Congo, Dem. Rep."	
			replace country = "Congo"	    						if country=="Congo, Rep."
			replace country = "Côte d'Ivoire"	    					if country=="Cote d'Ivoire"
			replace country = "Yemen"	    						if country=="Yemen, Rep."
			replace country = "Saint Lucia"    						if country=="St. Lucia"		
			replace country = "Saint Vincent and the Grenadines"    			if country=="St. Vincent and the Grenadines"	
			replace country = "Sint Maarten"    						if country=="Sint Maarten (Dutch part)"	
			replace country = "Saint Kitts and Nevis"    					if country=="St. Kitts and Nevis"	
			replace country = "Saint Lucia"    						if country=="St. Lucia"
		  
		*Add in map data	
		mmerge country using `world'
			drop if _merge==1
			drop _merge
		
		*Create same categories as used for region map
		gen gdpcat = .
		replace gdpcat = 1  if gdp <2040.4  & gdp!=.
		replace gdpcat = 2  if gdp >2040.4  & gdp<3450.6
		replace gdpcat = 3  if gdp >3450.6  & gdp<4713.9
		replace gdpcat = 4  if gdp >4713.9  & gdp<5908.5
		replace gdpcat = 5  if gdp >5908.5  & gdp<7709.9
		replace gdpcat = 6  if gdp >7709.9  & gdp<10856.6
		replace gdpcat = 7  if gdp >10856.6 & gdp<17621.5
		replace gdpcat = 8  if gdp >17621.5 & gdp<26176.4
		replace gdpcat = 9  if gdp >26176.4 & gdp<33862.7
		replace gdpcat = 10 if gdp >33862.7 & gdp!=.	
	
		spmap gdpcat using "world-c.dta", id(id)	///
		  fcolor( "220 0 20" "208 0 47" "185 0 70"  "162 0 93" "140 0 100" "120 0 120" "112 0 143" "86 0 169" "64 0 191" "47 0 208"  )  clm(unique) ///
		  legend(label(2 "< $2000") label(3 "$2000-3500") label(4 "$3500-4700") label(5 "$4700-5900") label(6 "$5900-7700") label(7 "$7700-10,800")  ///
		  label(8 "$10,800-17,600")  label(9 "$17,600-26,200") label(10 "$26,200-33,900") label(11 "> $33,900") ///
		  label(12 "No data")  size(vsmall) title("GDP per capita", size(vsmall)))  ndfcolor("192 192 192") ndocolor("192 192 192") ///
		  ocolor(none ...) saving(g1, replace) 
			   
*Add WDI data to region data, create graph with regions and country data where regional data is not available

		mmerge country region using `region'
		
		replace gdpcatr = gdpcat if gdpcatr==.
		replace gdpcatr = 11- gdpcatr
		duplicates drop

		spmap gdpcatr using "world-c.dta", id(id)	///
		  fcolor( "47 0 208"  "64 0 191" "86 0 169" "112 0 143"  "120 0 120" "140 0 100" "162 0 93" "185 0 70"  "208 0 47" "220 0 20" )  clm(unique) ///
		  legend(label(11 "< $2000") label(10 "$2000-3500") label(9 "$3500-4700") label(8 "$4700-5900") label(7 "$5900-7700") label(6 "$7700-10,800")  ///
		  label(5 "$10,800-17,600")  label(4 "$17,600-26,200") label(3 "$26,200-33,900") label(2 "> $33,900") ///
		  label(12 "No data") size(vsmall) title("Regional GDP per capita", size(vsmall)))  ///
		  ocolor(none ...)   ndfcolor("192 192 192") ndocolor("192 192 192") saving(g1, replace)