{"cells":[{"metadata":{},"cell_type":"markdown","source":"# This notebook cleans up the cities data, selects only US cities, and finds lat/lng/county for each city. \n\nNote: the kernel won't run without setting BING_API_KEY to your key. "},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"!pip install geocoder\n!pip install uszipcode","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","collapsed":true,"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"import numpy as np\nimport matplotlib.pyplot as plt\nimport pandas as pd\nfrom scipy import stats\nfrom sklearn.preprocessing import RobustScaler, QuantileTransformer\nimport geopandas as gpd\nimport geocoder\nimport time \nfrom tqdm import *\nimport plotly.express as px\n\nfrom uszipcode import Zipcode\nfrom uszipcode import SearchEngine\nsearch = SearchEngine(simple_zipcode=False)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let's open the cities disclosing file and select the cities in the US. "},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_disc_2020 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Disclosing/2020_Cities_Disclosing_to_CDP.csv\")\ncities_disc_2020_sub = cities_disc_2020[cities_disc_2020[\"Country\"] == \"United States of America\"]\n\norgs_2020 = list(cities_disc_2020_sub[\"Organization\"]) # list of organization names in US\ncities_2020 = list(cities_disc_2020_sub[\"City\"]) # list of cities in US\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"BING_API_KEY = \"YOUR_API_KEY\"\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let's do a Bing geocoding query for all the cities. There were a few quirks but the following conditionals seem to work! "},{"metadata":{"trusted":true},"cell_type":"code","source":"gs_2020 = []\nfor i in tqdm(range(len(orgs_2020))):\n    \n    if(orgs_2020[i] == \"City of Toledo\"): #toledo and birmingham need the \"city of\" in their string to return the right result\n        g = geocoder.bing(orgs_2020[i], key=BING_API_KEY)\n        gs_2020.append(g.raw)\n    elif(orgs_2020[i] == \"City of Birmingham\"):\n        g = geocoder.bing(orgs_2020[i], key=BING_API_KEY)\n        gs_2020.append(g.raw)\n    elif(orgs_2020[i][0:2] == \"Ci\"): #strip city of \n        g = geocoder.bing(orgs_2020[i][8:], key=BING_API_KEY)\n        gs_2020.append(g.raw)\n    else:\n        g = geocoder.bing(orgs_2020[i], key=BING_API_KEY)\n        gs_2020.append(g.raw)\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Okay, now we want to find the counties that span a city. Sometimes a city can have multiple counties (ex. NYC), so we will also be querying the population of city in each county and weight the corresponding county accordingly. The code below does that. It's a little messy but works. "},{"metadata":{"trusted":true},"cell_type":"code","source":"counties_2020 = [] # all the counties for a city\ncounties_weights_2020 = [] # \"weight\" each county by the population of the city that is in the county. this is done using uszipcode.\nbbox_2020 = [] # bbox of city. unreliable? \nloc_2020 = [] # latlng\nstates_2020 = [] # state \ncities_bing_name_2020 = [] # bing name more reliable\npop_densities_2020 = [] # pop_density\nmatch_multiple_county = np.zeros(len(cities_2020), dtype=bool)\nzip_2020 = []\n\nfor i in tqdm(range(len(cities_2020))):\n    split_city = ((gs_2020[i][\"address\"][\"formattedAddress\"].split(\", \")))\n\n    if(orgs_2020[i] == \"New York City\"): # NYC needs to be handcoded; data using Wiki\n        counties_2020.append('[\"New York County\", \"Kings County\", \"Bronx County\", \"Richmond County\", \"Queens County\"]')\n        counties_weights_2020.append([1.628/8.336, 2.559/8.336,1.412/8.336,0.476/8.336, 2.253/8.336])\n        result = search.by_city_and_state(split_city[0], split_city[1], returns=-1)\n        pop_densities_2020.append(result[0].population_density)\n        match_multiple_county[i] = True\n        zip_2020.append(result[0].zipcode)\n        \n    elif(orgs_2020[i] == \"District of Columbia\"): # DC also needs to be handcoded\n        counties_2020.append(orgs_2020[i])\n        counties_weights_2020.append(1.0)\n        pop_densities_2020.append(np.nan)\n        zip_2020.append(20001)\n        \n    else:\n        if(len(split_city) == 2): # it is a city! because the split led to something like [\"Seattle\", \"WA\"]\n\n            result = search.by_city_and_state(split_city[0], split_city[1], returns=-1) # search all zipcodes belonging to city\n            \n\n            if(len(result)==0): # if this does not work use Bing \n\n                counties_2020.append(gs_2020[i][\"address\"][\"adminDistrict2\"])\n                counties_weights_2020.append(1.0)\n                pop_densities_2020.append(np.nan)\n                zip_2020.append(np.nan)\n\n            else: # if this does work, we now have a list of zipcodes and their respective county belonging to a city\n                zip_2020.append(result[0].zipcode)\n                \n                all_counties = [] # get all the counties\n                for j in range(len(result)):\n                    if(len(result[j].county)>0):\n                        all_counties.append(result[j].county)\n    \n                all_counties = np.array(all_counties)\n                unique_counties, unique_counts = np.unique(all_counties, return_counts=True) # get unique counties bc zipcodes can belong to same county\n            \n                if(len(unique_counties) == 1): # only one county for this city\n                    counties_2020.append(unique_counties[0])\n                    counties_weights_2020.append(1.0)\n\n                    if(result[0].population_density == None):\n                        pop_densities_2020.append(np.nan)\n                    else:\n                        pop_densities_2020.append(result[0].population_density)\n\n                else:\n\n                    pops = np.array([result[j].population if result[j].population != None else 0 for j in range(len(result))]) # get population for this zip code\n                    pops_county = np.array([np.sum(pops[all_counties == unique_counties[j]]) for j in range(len(unique_counties))]) # combine population by county\n                    pop_density = np.array([result[j].population_density if result[j].population_density != None else 0 for j in range(len(result))]) # population density\n                    pop_density_county = np.array([np.mean(pop_density[all_counties == unique_counties[j]]) for j in range(len(unique_counties))]) # population density by ocunty\n                    weights = (pops_county/(pops_county.sum())) # county weights by population\n                    pop_density_city = np.sum(pop_density_county*weights)\n                    \n                    \n                    counties_2020.append(unique_counties.tolist())\n                    counties_weights_2020.append(weights.tolist())\n                    pop_densities_2020.append(pop_density_city)\n                    match_multiple_county[i] = True\n        else: # use bing\n            counties_2020.append(gs_2020[i][\"address\"][\"adminDistrict2\"])\n            counties_weights_2020.append(gs_2020[i][\"address\"][\"adminDistrict2\"])\n            pop_densities_2020.append(np.nan)\n            zip_2020.append(np.nan)\n    \n    bbox_2020.append(gs_2020[i][\"bbox\"])\n    loc_2020.append(\"POINT (%f %f)\"%(gs_2020[i][\"point\"][\"coordinates\"][1], gs_2020[i][\"point\"][\"coordinates\"][0]))\n    states_2020.append(gs_2020[i]['address'][\"adminDistrict\"])\n    cities_bing_name_2020.append(gs_2020[i][\"address\"][\"formattedAddress\"])\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Most datasets at the county level use the FIPS code as a unique identifier since a lot of counties can have the same name. We use the CDC SVI data to match the counties we derived above to their respective FIPS code"},{"metadata":{"trusted":true},"cell_type":"code","source":"df = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions//Supplementary Data/CDC Social Vulnerability Index 2018/SVI2018_US_COUNTY.csv\")\ncounties = list(cities_disc_2020_sub['counties'])\nstates = list(cities_disc_2020_sub['state'])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import ast \nmatch = []\n\nfor i in range(len(counties)): # loop over all derived counties\n    if(not match_multiple_county[i]): # only one county\n        if \"County\" in counties[i]:\n            counties[i] = counties[i].replace(\" County\",\"\")\n            counties[i] = ' '.join(counties[i].split())\n        elif \"Borough\" in counties[i]:\n            counties[i] = counties[i].replace(\" Borough\",\"\")\n            counties[i] = ' '.join(counties[i].split())\n        elif \"Parish\" in counties[i]:\n            counties[i] = counties[i].replace(\" Parish\",\"\")\n            counties[i] = ' '.join(counties[i].split())\n        else:\n            pass    \n    \n        idxs = np.where(df[\"COUNTY\"] == counties[i])[0]\n    \n            \n        if(len(idxs) == 1): # simple match \n            match.append(idxs[0])\n            \n        elif('DC' in counties[i]): # if it's DC \n            idxs = np.where(np.array([df[\"ST_ABBR\"] == \"DC\"]).flatten())[0]\n            match.append(idxs[0])\n            \n        else: # multiple counties share the same name, so we also gotta match it to the state\n            sts = df[\"ST_ABBR\"].iloc[idxs]\n            idxss = np.where(sts == states[i])[0]\n            match.append(idxs[idxss][0])\n    else: # has multiple counties\n        tmatch = []\n        mcounties = counties[i]\n        for j in range(len(mcounties)):\n            \n            if \"County\" in mcounties[j]:\n                mcounties[j] = mcounties[j].replace(\" County\",\"\")\n                mcounties[j] = ' '.join(mcounties[j].split())\n            elif \"Borough\" in mcounties[j]:\n                mcounties[j] = mcounties[j].replace(\" Borough\",\"\")\n                mcounties[j] = ' '.join(mcounties[j].split())\n            elif(mcounties[j] == \"Municipality of Anchorage\"):\n                mcounties[j] = mcounties[j].replace(\"Municipality of \",\"\")\n                mcounties[j] = ' '.join(mcounties[j].split())\n            elif \"Municipality\" in mcounties[j]:\n                mcounties[j] = mcounties[j].replace(\" Municipality\",\"\")\n                mcounties[j] = ' '.join(mcounties[j].split())                \n            elif \"Parish\" in mcounties[j]:\n                mcounties[j] = mcounties[j].replace(\" Parish\",\"\")\n                mcounties[j] = ' '.join(mcounties[j].split())\n            elif \"city\" in mcounties[j]:\n                mcounties[j] = mcounties[j].replace(\" city\",\"\")\n                mcounties[j] = ' '.join(mcounties[j].split())\n            else:\n                pass    \n\n            idxs = np.where(df[\"COUNTY\"] == mcounties[j])[0]\n            \n\n\n            if(len(idxs) == 1):\n                tmatch.append(idxs[0])\n            elif('DC' in mcounties[j]):\n                idxs = np.where(np.array([df[\"ST_ABBR\"] == \"DC\"]).flatten())[0]\n                tmatch.append(idxs[0])\n            elif(len(idxs) > 1):\n                sts = df[\"ST_ABBR\"].iloc[idxs]\n                idxss = np.where(sts == states[i])[0]\n                tmatch.append(idxs[idxss][0])\n            else:\n                tmatch.append(np.nan)\n        match.append(tmatch)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"match_fips = []\nflag_multiple = np.zeros(len(match), dtype=np.bool)\nfor i in range(len(match)):\n    if(type(counties[i]) == str or type(counties[i]) == np.str_):\n        match_fips.append(df[\"FIPS\"].iloc[match[i]])\n    else:\n        tm = []\n        for j in range(len(match[i])):\n            tm.append(df[\"FIPS\"].iloc[match[i][j]])\n        match_fips.append(tm)\n        flag_multiple[i] = True","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"a couple of the external datasets are also at the MSA level instead of county-level. We will choose a random zip code and derive the MSA name and MSA code using data from the US dept of labor (fs11_gpci_by_msa-ZIP.xls)"},{"metadata":{"trusted":true},"cell_type":"code","source":"zips_2020_unassigned = {\n    \"City of Piedmont, CA\" : \"94618\",\n    \"Orange County, NC\" : \"27243\",\n    \"Cuyahoga County\" : \"44118\",\n    \"Santa Fe County\" : \"87507\",\n    \"Boulder County\" : \"80302\",\n    \"Summit County, UT\" : \"84061\",\n    \"Dane County\" : \"53508\",\n    \"Broward County, FL\" : \"33024\",\n    \"City of Milwaukie, OR\" : \"97206\", \n    \"City and County of Honolulu\": \"96795\"\n}\n\nfor i in range(len(zip_2020)):\n    if(np.isnan(float(zip_2020[i]))):\n        zip_2020[i] = zips_2020_unassigned[orgs_2020[i]]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"zz = pd.read_excel(\"/kaggle/input/msabyzip/fs11_gpci_by_msa-ZIP.xls\", sheet_name=\"fs11gpci by ZIP, owcp\", skiprows=10)\n\nMSA = [] \nMSA_name = []\nfor i in range(len(zip_2020)):\n    cdf = zz[zz[\"ZIP CODE\"] == int(zip_2020[i])]\n    MSA.append(int(cdf[\"MSA No.\"]))\n    MSA_name.append(str(cdf[\"MSA Name\"].iloc[0]))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_disc_2020_sub[\"counties\"] = counties_2020\ncities_disc_2020_sub[\"counties_weights\"] = counties_weights_2020\ncities_disc_2020_sub[\"bbox\"] = bbox_2020\ncities_disc_2020_sub[\"loc\"] = loc_2020\ncities_disc_2020_sub[\"state\"] = states_2020\ncities_disc_2020_sub[\"city_bing\"] = cities_bing_name_2020\ncities_disc_2020_sub[\"pop_density\"] = pop_densities_2020\ncities_disc_2020_sub[\"FIPS\"] = match_fips\ncities_disc_2020_sub[\"flag_multiple_counties\"] = flag_multiple\ncities_disc_2020_sub[\"MSA\"] = MSA\ncities_disc_2020_sub[\"MSA_name\"] = MSA_name\ncities_disc_2020_sub[\"random_zipcode\"] = zip_2020","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_disc_2020_sub.to_csv(\"cities_updated_geo_us_2020.csv\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}