{"cells":[{"metadata":{},"cell_type":"markdown","source":"# **Tracking Tree Canopy Equality with Green Coverage Model**"},{"metadata":{},"cell_type":"markdown","source":"![](https://www.nationalgeographic.com/content/dam/science/2020/07/06/crown-shyness/crownshyness_mm9404_200619_00101.ngsversion.1594031243744.adapt.1900.1.jpg)\n*image source*: https://www.nationalgeographic.com/science/2020/07/tree-crown-shyness-forest-canopy/"},{"metadata":{},"cell_type":"markdown","source":"There are a multitude of benefits to tree coverage. Trees regulate extreme tempratures and manages flash flooding. Trees also adsorb carbon dioxide from the atmosphere, helping to offset greenhouse gas emissions. \n\nGiven these benefits, it is clear that trees have an implicit value. In this notebook we explore tree coverage inequaility within cities and how that perpetuates existing imbalances between the rich and the poor. The imbalances we discuss are public health outcomes related to extreme heat and higher vulnerability to flooding. No other category of hazardous weather event in the United States has caused more fatalities over the last few decades than extreme heat [1]\n\nThe KPI we suggest to CDP is the **Gini Coefficient** for a 20 km x 20 km region around a central coordinate. We calcuate this by querying satelitel images for the area of interest, then extract a subset of the green spectrum to find the density of green for the area. We then look at how evenly the green coverage is spread to calculate the distribution of green space. "},{"metadata":{},"cell_type":"markdown","source":"\n# Executive Summary\n\n<a id='1.1'></a>\n\n* Have created an automated tool which can take a city name as input to output a green coverage map and Lorenz Curve\n\n* Can provide actionable insights for any city by showing which areas are under represented with tree coverage and whether the tree distribution is equitable\n\n* Creates demand for low skilled work. A group sereverly impacted by the Covid-19 pandemic\n\n* Is an effective infrastructure project that can be started with relative ease and kick start an economy.  \n\n* Removes systematic bias towards minorities and lower income groups\n\n* As global temperatures and extreme weather rise, improves public health outcomes for low socio economic groups\n\n\nThe kernel is organised as follows. We start from the data description and proceed with the exploratory analysis. Next, a detailed method explanation is provided alongside examples of the application being used in Los Angeles and Sydney to ilustrate the robustness of the product. We then elaborate on the advantages and shortcomings of the approach and discuss further iterations of the tool."},{"metadata":{},"cell_type":"markdown","source":"# Exploratory Data Analysis on CDP Disclosure Data\n\nTo understand what KPIs would add value to cities we explored the disclosure data provided by CDP."},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"# Library Imports\n\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom IPython.core.display import HTML\nfrom IPython.display import  Markdown\nimport seaborn as sns\nimport random\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\nimport geopandas as gpd\nfrom shapely.geometry import Point\nfrom  sklearn.metrics import davies_bouldin_score\nfrom sklearn.decomposition import PCA\nfrom matplotlib.path import Path\nfrom matplotlib.spines import Spine\nfrom matplotlib.projections.polar import PolarAxes\nfrom matplotlib.projections import register_projection\nfrom matplotlib.patches import Circle, RegularPolygon\nfrom matplotlib.path import Path\nfrom matplotlib.projections.polar import PolarAxes\nfrom matplotlib.projections import register_projection\nfrom matplotlib.spines import Spine\nfrom matplotlib.transforms import Affine2D\nfrom sklearn import tree\nimport statsmodels.api as sm\nfrom scipy import stats\nimport json\nfrom sklearn.feature_extraction.text import CountVectorizer\nfrom sklearn.feature_extraction.text import TfidfVectorizer\nimport plotly.graph_objects as go\nimport json\nimport torch\n\n#centering figures\nHTML(\"\"\"<style> .output_png { display: table-cell; text-align: center; vertical-align: middle;}</style>\"\"\")\n\n#my colors\ncolors= ['#458B00','#629632','#397D02','#567E3A','#A6D785','#687E5A','#8AA37B','#476A34','#7BBF6A','#3D8B37','#426F42','#215E21']\n\n#remove pandas display limit\npd.options.display.max_colwidth = None","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"# Data Imports\n\ncities_2018 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2018_Full_Cities_Dataset.csv\",\n                         usecols=[c for c in list(pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2018_Full_Cities_Dataset.csv\", nrows =1)) if c not in ['Questionnaire','Year Reported to CDP','File Name','Last update','Comments']])\ncities_2019 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2019_Full_Cities_Dataset.csv\",\n                         usecols=[c for c in list(pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2019_Full_Cities_Dataset.csv\", nrows =1)) if c not in ['Questionnaire','Year Reported to CDP','File Name','Last update','Comments']])\ncities_2020 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2020_Full_Cities_Dataset.csv\",\n                         usecols=[c for c in list(pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2020_Full_Cities_Dataset.csv\", nrows =1)) if c not in ['Questionnaire','Year Reported to CDP','File Name','Last update','Comments']])\n\ncompanies_climate_change_2018 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2018_Full_Climate_Change_Dataset.csv\",\n                         usecols=[c for c in list(pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2018_Full_Climate_Change_Dataset.csv\", nrows =1)) if c not in ['Questionnaire','Year Reported to CDP','File Name','Last update','Comments']])\ncompanies_climate_change_2019 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2019_Full_Climate_Change_Dataset.csv\",\n                         usecols=[c for c in list(pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2019_Full_Climate_Change_Dataset.csv\", nrows =1)) if c not in ['Questionnaire','Year Reported to CDP','File Name','Last update','Comments']])\ncompanies_climate_change_2020 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2020_Full_Climate_Change_Dataset.csv\",\n                         usecols=[c for c in list(pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2020_Full_Climate_Change_Dataset.csv\", nrows =1)) if c not in ['Questionnaire','Year Reported to CDP','File Name','Last update','Comments']])\n\ngeo_cities_2020 = pd.read_csv(\"/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Disclosing/2020_Cities_Disclosing_to_CDP.csv\", usecols=[\"Account Number\",\"City Location\", \"Organization\", \"City\"]) ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"First we analysed what the biggest climate hazards facing cities today. As shown below flooding, extreme heat and extreme precipitation account for over 50% of the hazards reported."},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Isolating Cities Question\n\nquestion_3_responses = cities_2020[\n    (cities_2020['Question Number'] == '2.1') &\n    (cities_2020['Column Name'] == 'Climate Hazards')\n]\n\nquestion_3_responses = question_3_responses['Response Answer'].str.split(\">\", n = 1, expand = True).loc[:,0]\nquestion_3_response_count = question_3_responses.value_counts()\n\n# Isolating Coorperations Question\n\nquestion_2_responses = companies_climate_change_2020[\n    (companies_climate_change_2020['question_number'] == 'C2.3a') &\n    (companies_climate_change_2020['column_number'] == 3)\n]\n\nrisk_types = [\n    'Acute physical',\n    'Emerging regulation',\n    'Current regulation',\n    'Chronic physical',\n    'Market',\n    'Reputation',\n    'Technology'\n]\n\nquestion_2_responses = question_2_responses[~question_2_responses['response_value'].isin(risk_types)]\nquestion_2_response_count = question_2_responses['response_value'].value_counts().head(10)\n\n# Graphing\n\nfig = plt.figure(figsize=(15,9))\nax1 = fig.add_subplot(121)\n\nquestion_3_response_count.plot.pie(\n    textprops={'color':\"w\"},\n    pctdistance=0.7,\n    autopct='%.2f%%',\n    colors=colors, \n    labels=None,\n    ax=ax1,\n    ylabel=\"\"\n)\n\n\nax1.title.set_text(\"Climate Hazards Identified By Cities\")\nax1.legend(\n    question_3_response_count.index,\n    loc=\"lower center\", \n    bbox_to_anchor=(0.5, -0.4)\n)\n\nplt.show()\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Our next step was to understand what actions and goals cities are implementing to increase resilience to these hazards. Vectorising responses to Question 3.3 (2020) \n\n“Please describe the main goals of your city’s adaptation efforts and the metrics / KPIs for each goal.” \n\ngave us the results you see below. Interestingly the two most common actions were  ‘tree planting creation green’ and ‘planting creation green space’ In order to understand why these actions are so popular we then evaluated which climate issues the discussed actions addressed:"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"\npd.set_option('mode.chained_assignment', None)\n\nadaptation_actions_ques_3 = cities_2020[\n    (cities_2020['Question Number'] == '3.0') &\n    (cities_2020['Column Number'] == 2)\n]\n\nadaptation_actions_ques_3.loc[:, ['Response Answer']] = adaptation_actions_ques_3.loc[:, ['Response Answer']].fillna('No Response')\n\ndef get_top_n_bigram(corpus, n=None):\n    vec = TfidfVectorizer(ngram_range=(4,4), stop_words='english').fit(corpus)\n    bag_of_words = vec.transform(corpus)\n    sum_words = bag_of_words.sum(axis=0) \n    words_freq = [(word, sum_words[0, idx]) for word, idx in vec.vocabulary_.items()]\n    words_freq =sorted(words_freq, key = lambda x: x[1], reverse=True)\n    return words_freq[:n]\ncommon_words = get_top_n_bigram(adaptation_actions_ques_3['Response Answer'], 20)\n\nadaptation_actions_ques_3 = pd.DataFrame(common_words, columns = ['word' , 'count'])\n  \nadaptation_actions_ques_3  = adaptation_actions_ques_3.set_index('word')[:7]  ## Taking the first 5\n\n\n\nadaptation_goals_ques_3 = cities_2020[\n    (cities_2020['Question Number'] == '3.3') &\n    (cities_2020['Column Number'] == 1) \n]\nadaptation_goals_ques_3.loc[:, ['Response Answer']] = adaptation_goals_ques_3.loc[:, ['Response Answer']].fillna('No Response')\n\ndef get_top_n_bigram(corpus, n=None):\n    vec = TfidfVectorizer(ngram_range=(4,4), stop_words='english').fit(corpus)\n    bag_of_words = vec.transform(corpus)\n    sum_words = bag_of_words.sum(axis=0) \n    words_freq = [(word, sum_words[0, idx]) for word, idx in vec.vocabulary_.items()]\n    words_freq =sorted(words_freq, key = lambda x: x[1], reverse=True)\n    return words_freq[:n]\ncommon_words = get_top_n_bigram(adaptation_goals_ques_3['Response Answer'], 20)\n# for word, freq in common_words:\n#     print(word, freq)\nadaptation_goals_ques_3 = pd.DataFrame(common_words, columns = ['word' , 'count'])\n\nadaptation_goals_ques_3.loc[adaptation_goals_ques_3.word=='reforzar el sistema salud','word']= 'make a stronger health care'\nadaptation_goals_ques_3.loc[adaptation_goals_ques_3.word=='monitoramento risco em tempo','word']= 'monitoring risk in real time'\n\nadaptation_goals_ques_3  = adaptation_goals_ques_3.set_index('word')[:6]  ## Taking the first 5\n\n# Graphing\n\nfig = plt.figure(figsize=(15,9))\nax1 = fig.add_subplot(121)\nax2 = fig.add_subplot(122)\n\nadaptation_actions_ques_3['count'].plot.pie(\n    textprops={'color':\"w\"},\n    pctdistance=0.7,\n    explode=(0.1,0.1,0,0,0,0,0),\n    autopct='%.2f%%',\n    colors=colors, \n    labels=None,\n    ax=ax1,\n    ylabel=\"\"\n)\n\nadaptation_goals_ques_3['count'].plot.pie(\n    textprops={'color':\"w\"},\n    explode=(0.1,0,0,0,0,0),\n    pctdistance=0.7,\n    autopct='%.2f%%',\n    colors=colors, \n    labels=None,\n    ax=ax2,\n    ylabel=\"\"\n)\n\nax1.title.set_text(\"Actions - Present\") ## present \nax2.title.set_text(\"Goals/KPI's - Future\") ## future\nax1.legend(\n    adaptation_actions_ques_3.index,\n    loc=\"lower center\", \n    bbox_to_anchor=(0.5, -0.4)\n)\nax2.legend(\n    adaptation_goals_ques_3.index,\n    loc=\"lower center\", \n    bbox_to_anchor=(0.5, -0.4)\n)\n\nplt.show()\n\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"As the above visualization shows, improving green coverage increases resilience towards each of the top three risks we discovered earlier, clearly justifying its popularity."},{"metadata":{},"cell_type":"markdown","source":"By selecting these actions we were then able to observe the metrics used to measure their implementation and success. As shown below, this analysis showed us that cities do not have a common KPI to measure the same objective. \n\nUltimately our processing and visualisation of CDP data revealed the demand for a green coverage KPI that could be easily used by any city, anywhere in the world. Such a tool would allow cities around the world to compare their urban greenery and its social distribution to other cities. As a global organisation who has set themselves apart in disclosure we believe this KPI is perfect for CDP."},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"question_3_responses = cities_2020[cities_2020['Question Number'] == '3.3'].copy()\n\nquestion_3_3_responses = cities_2020[(cities_2020['Question Number'] == '3.3') &\n                                     (cities_2020['Column Number'] == 1)].copy()\n\nquestion_3_3_responses['Response Answer'] = question_3_3_responses['Response Answer'].fillna('No Response')\n \n# Filter question 3.3 to cities that mention green spaces in adaptation goals\n\n# list_of_strings = ['tree', 'plant', 'park', 'planting', 'canopy', 'canopies']\nlist_of_strings = ['tree']\n\n\npattern = '|'.join(list_of_strings)\n\nquestion_3_3_filtered = question_3_3_responses[question_3_3_responses['Response Answer'].str.contains(pattern)]\n\n# Retrieve all columns for the rows that have the list of strings in thier adaptation goals.\n\nall_green_responses = pd.DataFrame()\n\nfor index, row in question_3_3_filtered.iterrows():\n\n    all_green_responses = all_green_responses.append(question_3_responses[(question_3_responses['Account Number'] == row['Account Number']) &\n                                          (question_3_responses['Row Number'] == row['Row Number'])], ignore_index=True)\n    \nall_green_responses_2 = all_green_responses['Response Answer'][all_green_responses['Column Number'] == 2]\n\nall_green_responses_2 = all_green_responses_2.str.split(\">\", n=1, expand=True).loc[:, 0].value_counts()\n\nfig = go.Figure(data=[go.Sankey(\n    node = dict(\n      pad = 15,\n      thickness = 20,\n      line = dict(color = \"#215E21\", width = 0.5),\n      label = [\"Tree\", \"Extreme hot temperature\", \"Flood and sea level rise\", \"Extreme Precipitation\", \"Storm and wind\", \"Water Scarcity\", \"Chemical change\", \"Extreme cold temperature\", \"Biological hazards\", \"Mass movement\"],\n      color = \"#4BB74C\"\n    ),\n    link = dict(\n      source = [0, 0, 0, 0, 0, 0, 0, 0, 0], # indices correspond to labels, eg A1, A2, A1, B1, ...\n      target = [1, 2, 3, 4, 5, 6, 7, 8, 9],\n      value = [20, 10, 8, 6, 5, 3, 2, 2, 1],\n    color = \"#C1FFC1\"\n  ))])\n\nfig.update_layout(title_text=\"Hazard related to Tree\", font_size=15)\nfig.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"\n# The intersection between the tree coverage and social inequality\n\n## Health outcomes \n\nMedical and health researchers have shown that fatalities during heat waves are most commonly due to respiratory and cardiovascular diseases, primarily from heat's negative effect on the cardiovascular system. In an attempt to control one's internal temperature, the body’s natural instinct is to circulate large quantities of blood to the skin. However, to perform this protective measure against overheating actually harms the body by inducing extra strain on the heart. This excess strain has the potential to trigger a cardiac event in those with chronic health problems, such as the elderly, (Cui et al.) Frumkin showed that the relationship of mortality and temperature creates a J-shaped function, showing a steeper slope at higher temperatures. Records show that more casualties have resulted from heat waves than hurricanes, floods, and tornadoes together. This statistic’s significance is that extreme heat events (EHEs) are becoming more frequent (Stone et al)\n\nStudies held from 1989 to 2000 have also recorded a rise of 5.7% in mortality during heat waves. A study revealed that Rome’s elderly population endures a higher mortality rate during heat waves, at 8% excess for the 65–74 age group and 15% for above 74 (Schifano et al.). Another study found that in French cities during the 2003 heat wave, small towns saw an average excess mortality rate of 40%, while Paris witnessed an increase of 141%. During this period, a 0.5 °C increase above the average minimum nighttime temperature doubled the risk of death in the elderly (Dousset et al.).\n\nSince the air temperature of urban areas with more trees can be around 4C cooler than those without. On a more local level, the air temperature on a treeless residential street can be 10C higher than a nearby shady street (Sinfield et al). \n\nA city’s poorest areas tend to have less tree canopy than wealthier areas, a pattern that is especially pronounced on the concrete-dense neighbourhoods where temperature regualtion is most important (Ready et al.).\n\nThe areas which have faced systematic inequality with tree coverage are now seeing poorer public health outcomes because of it. Through no fault of their own, at risk people are now more vulnerable to extreme heat because of where they live.\n\nThe best way to counter this issue is to plan trees in the areas which need it most. By planting trees, city councils can effectively regular temperatures and therefore improve public health outcomes for vulnerable memebers of thier society. Since global temperatures are rising, this is an issue which will only become more prevelant as time goes on. It is critical for councils to act now, which is what has inspired our KPI.\n\n\n## Danger from flooding\n\nThere are several types of flooding that affect communities:\n-         River Flooding: When the amount of water entering a river exceeds its holding capacity and overflows its banks\n-         Surface water flooding – when heavy rainfall runs on hard or saturated surfaces without getting absorbed into the ground and collects in low areas damaging properties and misplacing communities\n-         Drain and sewer flooding – When heavy rain causes drain and sewers to be blocked\n-         Coastal flooding – due to climate change, weather and tidal conditions cause an increase in sea levels affecting coastal neighbourhoods.\n\nPlanting trees and extending greenspaces helps by:\nCanopy cover – Rainfall is intercepted by trees significantly slowing its speed. This causes about 30% of the water to evaporate back into the atmosphere without reaching the ground. This can even occur in winter where trees intercept and re-evaporate rainfall\nRoots – The trees’ root systems allow water to penetrate deeper into the soil which decreases surface run-off while also increasing water storage capacity of the soil\nUrban trees – In urban spaces, an increase in impermeable surfaces such as roads, pavements and driveways has led to increased surface water run-off. Trees reduce surface run-off by 80% compared to asphalt. Incorporating green spaces in urban spaces would drastically reduce run-off leading to a decreased risk of flooding (Woodland). \n\n\n## How can this pull cities out of a recession without perpetuating socical inequities?\n\nTree planting projects can help pull cities out of a recession by moblising the hardest hit workers from the Covid-19 pandemic, low skilled workers. By getting this group employeed, the economy can pull itself out of the recession as it will put money into the hands of many people. This will stimulate local economies and get the gears moving again. Being a labour intesive project, with virtually no qualifications required, this is a great way to get the masses employeed in any city in the world. Also, trees often need at least two years of costant care (Ready et al.) which means these employees can stay hired until the economy has rebounded and the labour market is thriving again. Now contrast this with investment in renewable energy. Of course, renewables are great, but the workers involved are typically high skilled and are most likely still in the job market and if not, there is plenty of oppoturnitites out there for them. By investing in renewables we are only worsening the income inequality as the low skilled workers are being forgotten about. \n\n\n## How can this be done in a socialy equitable way with the backdrop of a global pandemic?\nAs mentioned above, by mobilising low skilled workers for these projects, cities are able to put money in the hands of people who need ti most. This contributes to making society more equitable since it is evens out the income distribution. With industries and hospitality practically reaching a stand still with the pandemic, we can help these newly unemployeed workers find meaningful work again.\n\n\n## How can corporations help solve this problem?\nCities could outsource the labour requirement to corporations as it is not a councils expertise to directly hire and manage a workforce like what is anticipated. \n\n\n## How does this measure the intersection in the context of resiliency?\nThe intersection allows cities to be resilient in a mryiad of different ways. By planting more trees, city's are defending themselves against rising temperatures, they are protecting their most vulnerable citizens and are extracting CO2 from the atmosphere as trees are able to absorb it.\n\n"},{"metadata":{},"cell_type":"markdown","source":"# Tree Canopy Equality KPI\n\nOur KPI allows cities around the world to evaluate their green coverege in a consistent and reliable way. We create an interface for collecting and assesing satelite images around the world and the equality of greenery distribution.\n\n## Method:\n\nWe follow the following steps:\n\n1. Get coordinates of city center from CDP data\n2. Fetch zoomed out satelite images\n3. Determine city boundries using custom ML model\n4. Collect thousands of zoomed in images within the city bounds and asses green coverage\n5. Combine results to form lorenz curve and calculate gini coefficient\n6. Visualise areas with least tree coverage on heat map\n\nCode:\n\nWe have encapsulated this logic in the following blocks and included them below:\n1. Geographical grid search system\n2. Machine learning based city detection\n3. Satelite imagery pipeline\n4. Complete pipeline\n"},{"metadata":{},"cell_type":"markdown","source":"## Green Coverage Tool"},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"# 1. Geographical grid search system \n\ndef getBoundingBoxes(bounds, zoom=14):\n    \n    resolution_map = {\n        14: 0.054885,\n        18: 0.003251\n    }\n    resolution = resolution_map[zoom]\n\n    xSorted = sorted(bounds, key=lambda p: p[0])\n    minX = xSorted[0][0]\n    maxX = xSorted[-1][0]\n\n    ySorted = sorted(bounds, key=lambda p: p[1])\n    minY = ySorted[0][1]\n    maxY = ySorted[-1][1]\n\n\n    gridWidth = int((maxX - minX) / resolution)\n    gridHeight = int((maxY - minY) / resolution)\n    return [\n        [\n            [minX + x * resolution, minY + y * resolution],\n            [minX + (x + 1) * resolution, minY + y * resolution],\n            [minX + (x + 1) * resolution, minY + (y + 1) * resolution],\n            [minX + x * resolution, minY + (y + 1) * resolution],\n            [minX + x * resolution, minY + y * resolution],\n        ]\n        for x in range(gridWidth)\n        for y in range(gridHeight)\n    ]","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"# 2. Machine learning based city detection\n# Only the functions needed for inference are included within this document\n# All training was completed in a seperate directory to avoid complexity balooning\n\nclass BasicClassificationDataset(torch.utils.data.Dataset):\n    def __init__(self, dir=None, image_paths=None):\n        self.dir = dir\n        self.image_paths = image_paths\n        self.transforms = transforms.Compose(\n            [\n                transforms.ToTensor(),\n                transforms.Normalize((0.5, 0.5, 0.5), (0.5, 0.5, 0.5)),\n                # transforms.RandomHorizontalFlip(0.5),\n                # transforms.RandomVerticalFlip(0.5),\n            ]\n        )\n        if dir:\n            self.labels = os.listdir(dir)\n            self.images = []\n\n            for i in range(len(self.labels)):\n                label = self.labels[i]\n                images = os.listdir(os.path.join(dir, label))\n                self.images += [\n                    [os.path.join(dir, label, image), i] for image in images\n                ]\n        if image_paths:\n            self.labels = [\"None\"]\n            self.images = list(map(lambda p: [p, 0], self.image_paths))\n\n    def __getitem__(self, idx):\n        image_path, image_class = self.images[idx]\n        img = Image.open(image_path).convert(\"RGB\")\n        img = self.transforms(img)\n        return img, image_class\n\n    def __len__(self):\n        return len(self.images)\n        \ndef urban_classification(image_paths):\n    dataset = BasicClassificationDataset(image_paths=image_paths)\n    loader = torch.utils.data.DataLoader(\n        dataset, batch_size=1, shuffle=False, num_workers=0\n    )\n\n    classifications = []\n\n    with torch.no_grad():\n        for i, data in enumerate(loader, 0):\n            images, labels = data\n            images = images.to(device)\n            outputs = model(images)\n            _, predicted = torch.max(outputs.data, 1)\n            classifications += predicted.tolist()\n\n    return classifications\n\n## Basic color detection for trees\n## If given more time we would train a segmentation model for the task\n\ndef get_green_mask(image_path):\n    img = cv2.imread(image_path)\n    img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n    width, height = img.shape[:-1]\n\n    mask = np.zeros([width, height])\n\n    for x in range(width):\n        for y in range(height):\n            pixel = img[x][y]\n            if (pixel[1]-2) > pixel[0] and (pixel[1] - 2) > pixel[2] and (pixel[1]<200):\n                mask[x][y] = 1\n\n    green_score = mask.sum() / (width * height)\n\n    return (img, mask, green_score)\n\ndef display_mask(image_path):\n    img, mask, green_score = get_green_mask(image_path)\n\n    print(green_score)\n    fig = plt.figure(figsize=(15, 9))\n    ax1 = fig.add_subplot(121)\n    ax2 = fig.add_subplot(122)\n\n    ax1.imshow(mask)\n    ax2.imshow(img)\n    plt.show()\n    ","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"# 3. Satelite imagery pipeline\n\ndef satelite_download(tile, output_dir, zoom):\n    size = [600, 600]\n    zoom = zoom\n\n    centerX = (tile[0][0] + tile[2][0]) / 2\n    centerY = (tile[0][1] + tile[2][1]) / 2\n\n    image_url = f\"http://maps.googleapis.com/maps/api/staticmap?center={centerY},{centerX}&zoom={zoom}&size={size[0]}x{size[1]}&maptype=satellite&key={API_key}\"\n    filename = f\"{output_dir}/{centerX},{centerY}.png\"\n\n    if os.path.exists(filename):\n        return filename\n    else:\n        # Open the url image, set stream to True, this will return the stream content.\n        r = requests.get(image_url, stream=True)\n\n        # Check if the image was retrieved successfully\n        if r.status_code == 200:\n            # Set decode_content value to True, otherwise the downloaded image file's size will be zero.\n            r.raw.decode_content = True\n\n            # Open a local file with wb ( write binary ) permission.\n            with open(filename, \"wb\") as f:\n                shutil.copyfileobj(r.raw, f)\n\n            print(\"Image sucessfully Downloaded: \", filename)\n            return filename\n        else:\n            print(\"Image Couldn't be retreived\")\n            return None","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"# 4. Complete pipeline\n\ndef get_green_scores(city_name, API_key, data_dir, output_dir):\n\n    if (API_key=='') or (API_key==None):\n        print('Google Maps API Key Needed')\n    \n    device = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n    MODEL_PATH = \"drive/My Drive/Colab Notebooks/pleasework.pth\"\n    model=torch.load(MODEL_PATH)\n    model.eval()\n\n    longitude, latitude  = geo_cities_2020[geo_cities_2020[\"City\"] == \"Los Angeles\"]\\\n        [\"City Location\"].iloc[0].replace('POINT (', '').replace(')', '').split(\" \")\n    longitude, latitude = float(longitude), float(latitude)\n\n    width = 0.7\n    city_surrounding_area = [\n        [longitude - width/2,latitude - width/2],\n        [longitude + width/2,latitude + width/2]\n    ]\n    large_grids = getBoundingBoxes(city_surrounding_area, zoom=14)\n    large_image_paths = [download(grid, data_dir, 14) for grid in large_grids]\n\n    urban_grid_mask = urban_classification(large_image_paths)\n    urban_grids = [grids[i] for i in range(len(urban_grid_mask)) if (urban_grid_mask[i]==0)]\n\n    small_grids = []\n    for grid in large_grids:\n        small_grids += getBoundingBoxes(grid, zoom=18)\n    small_image_paths = [download(grid, data_dir, 18) for grid in small_grids]\n\n    for path in small_image_paths:\n        try:\n            file_name = path.split(\"/\")[-1]\n        except:\n            print(\"Something went wrong\")\n\n        if (count % 10 == 0):\n            with open(f\"{output_dir}/green_scores.json\", \"w\") as outfile:\n                json.dump(mapped_data, outfile)\n                print(\"saved\")\n\n\n        if file_name not in mapped_data:\n            centerX = float(file_name.split(\",\")[1].replace(\".png\", \"\"))\n            centerY = float(file_name.split(\",\")[0])\n            _, mask, green_score = get_green_mask(path)\n            mapped_data[file_name] = green_score\n\n    with open(f\"{output_dir}/green_scores.json\", \"w\") as outfile:\n        json.dump(mapped_data, outfile)\n\n    return mapped_data","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Example: Los Angeles"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"city_name = 'Los Angeles'\nAPI_key = ''\ndata_dir = './data'\noutput_dir = './outputs'\n\n# Here we display data that has been calculated through the following command\n# get_green_scores(city_name, API_key, data_dir, output_dir)\n\nwith open(\"/kaggle/input/green-score-tree-canopy/green_scores_la.json\") as json_file:\n    la_mapped_data = json.load(json_file)\n\nlowest_green_scores = []\nfor file_name in la_mapped_data:\n    centerX = float(file_name.split(\",\")[1].replace(\".png\", \"\"))\n    centerY = float(file_name.split(\",\")[0])\n    if la_mapped_data[file_name]<0.005:\n        lowest_green_scores.append([centerX, centerY, la_mapped_data[file_name]])\n\nimport folium\nfrom folium.plugins import HeatMap\n\nbase_map = folium.Map(location=[34.03273949999999, -118.3654856], control_scale=True, zoom_start=12, max_zoom=13, min_zoom=10)\nHeatMap(data=lowest_green_scores, radius=10, max_zoom=20, gradient={0.0: 'red', 0.05: 'red'}).add_to(base_map)\nbase_map","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Example: Sydney"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"city_name = 'Sydney\nAPI_key = ''\ndata_dir = './data'\noutput_dir = './outputs'\n\n# Here we display data that has been calculated through the following command\n# get_green_scores(city_name, API_key, data_dir, output_dir)\n\nwith open(\"/kaggle/input/green-score-tree-canopy/green_scores_sydney.json\") as json_file:\n    syd_mapped_data = json.load(json_file)\n\n    \nlowest_green_scores = []\nfor file_name in syd_mapped_data:\n    centerX = float(file_name.split(\",\")[1].replace(\".png\", \"\"))\n    centerY = float(file_name.split(\",\")[0])\n    if syd_mapped_data[file_name]<0.1:\n        lowest_green_scores.append([centerX, centerY, syd_mapped_data[file_name]])\n\n        import folium\nfrom folium.plugins import HeatMap\n\nbase_map = folium.Map(location=[-33.9358211, 151.03409150000005], control_scale=True, zoom_start=12, max_zoom=13, min_zoom=10)\nHeatMap(data=lowest_green_scores, radius=10, max_zoom=20, gradient={0.0: 'red', 0.05: 'red'}).add_to(base_map)\nbase_map","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Lorenz Curve and Gini Coefficient\n\nTo allow cities to measure their overall Tree Canopy Equality, we calculuate a Lorenz curve in realtime so councils can have get a graphical representation of what their unique distribution is like. This curve allows cities to see the proportion of overall tree canopy assumed by the bottom x% of the land area. \n*(Rememeber an individual land area is regarded at the 10 m x 10 m referenced above)*\n\nWe compute the Gini Coefficient by calculating the area between the line of equality and the lorenz curve. This final value is the **KPI** which cities want to optimise. \n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"with open(\"/kaggle/input/green-score-tree-canopy/green_scores_sydney.json\") as json_file:\n    syd_mapped_data = json.load(json_file)\nwith open(\"/kaggle/input/green-score-tree-canopy/green_scores_la.json\") as json_file:\n    la_mapped_data = json.load(json_file)\n\nsyd_ma = []\nla_ma = []\n\nfor file_name in syd_mapped_data:\n    centerX = float(file_name.split(\",\")[1].replace(\".png\", \"\"))\n    centerY = float(file_name.split(\",\")[0])\n    if syd_mapped_data[file_name]<10:\n        syd_ma.append(syd_mapped_data[file_name])\n\nfor file_name in la_mapped_data:\n    centerX = float(file_name.split(\",\")[1].replace(\".png\", \"\"))\n    centerY = float(file_name.split(\",\")[0])\n    if la_mapped_data[file_name]<10:\n        la_ma.append(la_mapped_data[file_name])\n  \n                 \n\ndef gini(arr):\n    ## first sort\n    sorted_arr = arr.copy()\n    sorted_arr.sort()\n    n = arr.size\n    coef_ = 2. / n\n    const_ = (n + 1.) / n\n    weighted_sum = sum([(i+1)*yi for i, yi in enumerate(sorted_arr)])\n    return coef_*weighted_sum/(sorted_arr.sum()) - const_\n\n\ndef lorenz_curve(arr, ax):\n    sorted_arr = arr.copy()\n    sorted_arr.sort()\n    X_lorenz = sorted_arr.cumsum() / sorted_arr.sum()\n    X_lorenz = np.insert(X_lorenz, 0, 0)\n    X_lorenz[0], X_lorenz[-1]\n    ## scatter plot of Lorenz curve\n    ax.scatter(np.arange(X_lorenz.size) / (X_lorenz.size - 1), X_lorenz,\n               marker='x', color='darkgreen', s=2)\n    ## line plot of equality\n    ax.plot([0, 1], [0, 1], color='k')\n    \nfig = plt.figure(figsize=(15,9))\nax1 = fig.add_subplot(121)\nax2 = fig.add_subplot(122)\n\n\nx = np.array(syd_ma)\n\nlorenz_curve(x, ax1)\n\ny = np.array(la_ma)\n\nlorenz_curve(y, ax2)\n\nax1.title.set_text(f\"Sydney - Gini Coefficient = {str(gini(x))[:5]}\") ## present \nax2.title.set_text(f\"Los Angeles - Gini Coefficient = {str(gini(y))[:5]}\") ## future","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Advantages\n\nThe best benchmark for our Green Coverage Tool and Gini Coefficent is Treepedia by MIT Senseable City Lab. They have tried to capture the green canopy by utilising Google Street view. They scan each panorama and use a predictive model to idennify tree canopy size from eye level. The outputs they provide is a map of the city where each street view location as a dot with a varying green intensity depending on how much tree was in the panorama (this shows a similar visualisation to our heat map). They also provide a Green View Index which is shows the percentage of canopy coverage for the entire city (similar to our gini coefficient). \n\nOur solution is able to recognise considerably more greenery, including park lands and backyards. As discussed these features increase ressilience to extreme weather events and therefore should be considered. Our product also has the advantage of directly integrating with CDPs data and customers."},{"metadata":{},"cell_type":"markdown","source":"## Disadvantages\n\nAlthough our work represents a great start there are still many improvements that can be made. These include:\n\n1. Training a ML segmentation model for identifying greenery rather than a color based approach\n2. Training another classification model to identify area type at a smaller scale\n3. Create an online interface that alows organisation reporting to CDP to easily view and interact with the data\n\nEach of these steps is possible and would greatly improve the value of our product. Time constraits prevented us from persuing them."},{"metadata":{},"cell_type":"markdown","source":"# Final Thoughts\n\nOverall this was a great exercise in predicting tree canopy inequality in cities. Although we see much more room for growth and improvement with our models and visualisations we belive that this can be a great tool that even now can be used by cities to see where to plant trees and to see how equal their current tree structure is.\n\nWe believe the generic nature of the tool strongly complements CDP's mission for global disclosure."},{"metadata":{},"cell_type":"markdown","source":"\n## Our team\n\nWe had the pleasure of working with an experienced and diverse team. Scattered around the world we connected over zoom, slack and git: \n\n<table>\n  <tr>\n    <td><img src=\"https://media-exp1.licdn.com/dms/image/C5603AQHbEQCGU89V6A/profile-displayphoto-shrink_400_400/0/1596957301931?e=1611792000&v=beta&t=fZn_O9amBy7bHGYLu9sd1XEIwjURey3TW-AWQ5Mtn1w\" width=\"150\"></td>\n      <td style='text-align: left;'><strong>Adrian Sarstedt</strong><br/><i>Data Engineer at the Florey Institute</i><br/><i>Australia</i></td>\n      <td><img src=\"https://media-exp1.licdn.com/dms/image/C5603AQHVTBlErFKspw/profile-displayphoto-shrink_400_400/0/1551742466431?e=1611792000&v=beta&t=wBO3P2BMWHC39dkXI3Kkjq5gKV4mZLoFeCZs20YaomQ\" width=\"150\"></td>\n    <td style='text-align: left;'><strong>Hamish Gunasekara</strong><br/><i>Data/Risk Analyst at Afterpay</i><br/><i>Australia</i></td>\n  </tr>\n    <tr>\n    <td><img src=\"https://media-exp1.licdn.com/dms/image/C4D03AQHeh45CPg4YsA/profile-displayphoto-shrink_400_400/0?e=1611792000&v=beta&t=IDlea0_mutyIpY6FKT4Btbtc3DtmAB8TFSoAkdfejWE\" width=\"150\"></td>\n    <td style='text-align: left;'><strong>Adham Al Hossary</strong><br/><i>Data Scientist at C-Capture</i><br/><i>England</i></td>\n      <td><img src=\"https://media-exp1.licdn.com/dms/image/C4D03AQFUX0Sfw87-JQ/profile-displayphoto-shrink_400_400/0?e=1611792000&v=beta&t=klWv7UqcQhpXRzP2tjgsLOWIdMYlFMd_ZHrCYyuutes\" width=\"150\"></td>\n    <td style='text-align: left;'><strong>Alexandra Golab</strong><br/><i>Business Analyst at CitiBank</i><br/><i>Spain</i></td>\n  </tr>\n  <tr>\n    <td><img src=\"https://media-exp1.licdn.com/dms/image/C4D35AQG0jqEpq15Lxw/profile-framedphoto-shrink_400_400/0?e=1606993200&v=beta&t=goN2AEOlVTO0ImAUhFJYXSzuvfzr5x3Ga-b-iobtyBQ\" width=\"150\"></td>\n    <td style='text-align: left;'><strong>Xavier Pivan</strong><br/><i>Deep Learning Engineer</i><br/><i>Portugal</i></td>\n  </tr>\n</table>"},{"metadata":{},"cell_type":"markdown","source":"# References\n\nThe study on heat and redlining:\nThe Effects of Historical Housing Policies on Resident Exposure to Intra-Urban Heat: A Study of 108 US Urban Areas - Jeremy S. Hoffman, Vivek Shandas and Nicholas Pendleton\nhttps://www.mdpi.com/2225-1154/8/1/12...\n\nInteractive maps of neighborhood heat and redlining:\nhttps://www.arcgis.com/apps/dashboard...\n\nRobert K. Nelson, LaDale Winling, Richard Marciano, Nathan Connolly, et al., “Mapping Inequality,” American Panorama, ed. Robert K. Nelson and Edward L. Ayers, accessed August 4, 2020.\nhttps://dsl.richmond.edu/panorama/red...\n\nAll about Urban Heat Islands from Climate Central [PDF]:\nhttp://assets.climatecentral.org/pdfs...\n\nCan trees and woods reduce flooding:\nhttps://www.woodlandtrust.org.uk/trees-woods-and-wildlife/british-trees/flooding/\n\nSatellite monitoring of summer heat waves in the Paris metropolitan area:\nhttps://rmets.onlinelibrary.wiley.com/doi/10.1002/joc.2222\n\nUrban tree canopy governance and redlined neighborhoods: an analysis of five cities\nhttps://dspace.mit.edu/handle/1721.1/127588\n\nResidential housing segregation and urban tree canopy in 37 US Cities\nhttps://osf.io/preprints/socarxiv/97zcs/\n"}],"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}