{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport folium\nimport seaborn as sns\nimport matplotlib as pl\n%matplotlib inline\nimport json\nimport os","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The objective of this challenge is to estimate an average historical emissions factor of electricity generated in a sub-national region - in this case, Puerto Rico.\n\nAn emission factor represents the amount of pollution released to the atmosphere in relation to the amount of the activity associated with those pollutants. In this case, the object of interest is electricity production. Thus, what we want to find out is how much pollution is released for each unit of electricity produced in Puerto Rico.\n\nThe formula according to the challenge documentation for total emissions is:\n\n> E = A x EF x (1-ER/100)\n\nwhere:\n\nE = emissions;\nA = activity rate;\nEF = emission factor, and\nER = overall emission reduction efficiency, %\n\nWhich means that:\n\n> EF = E / (A x (1-ER/100))"},{"metadata":{},"cell_type":"markdown","source":"Let's take a look at the data provided with the challenge. First, on the Global Power Plant database provided by WRI we can find data on all the power plants in Puerto Rico. From this dataset we can quickly compute A."},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"# Let's import plant data\ngp = pd.read_csv('/kaggle/input/ds4g-environmental-insights-explorer/eie_data/gppd/gppd_120_pr.csv')\n\n#Flatten geo column\ngp_f = gp.join(pd.io.json.json_normalize(gp[\".geo\"].map(json.loads).tolist()).add_prefix(\"geo.\")).drop([\".geo\"], axis=1)\ngp_f['longitude'] = gp_f['geo.coordinates'].map(lambda x: x[0])\ngp_f['latitude'] = gp_f['geo.coordinates'].map(lambda x: x[1])\ngp_f","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sw = [gp_f['latitude'].min(), gp_f['longitude'].min()]\nne = [gp_f['latitude'].max(), gp_f['longitude'].max()]\n\nm = folium.Map(tiles=\"Stamen Terrain\")\nm.fit_bounds([sw, ne])\n\nfor i in range(0, len(gp_f.index)):\n    folium.Marker([gp_f.iloc[i]['latitude'], gp_f.iloc[i]['longitude']], popup=gp_f.iloc[i]['primary_fuel']).add_to(m)\nm","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"As we can see, there are a total of 35 power plants in the region with a total of:"},{"metadata":{"trusted":true},"cell_type":"code","source":"print(\"Total capacity:\", gp_f['capacity_mw'].sum(), \"MW/h\")\nprint(\"Estimated cumulative generation:\", gp_f['estimated_generation_gwh'].sum(), \"GW/h\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plot = gp_f['primary_fuel'].value_counts().plot.pie(autopct='%1.0f%%')\nplot.set_title(\"Sources by # of plants\")\npl.pyplot.show()\n\nplot = gp_f['primary_fuel'].value_counts().plot.bar()\nplot.set_title(\"# of plants by source\")\npl.pyplot.show()\n\nprint(\"Top 3:\", gp_f[gp_f['primary_fuel'].isin(['Hydro', 'Gas', 'Solar'])]['capacity_mw'].sum() * 100 / gp_f['capacity_mw'].sum())\nprint(\"Renewables:\", gp_f[gp_f['primary_fuel'].isin(['Hydro', 'Wind', 'Solar'])]['capacity_mw'].sum() * 100/ gp_f['capacity_mw'].sum())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can see hydro, gas and solar dominating the number of plants, with 74% of plants coming from those sources, although they represent only 22% of the capacity. Oil, wind and coal follow. In total, renewables (hydro, solar and wind) amount to 6% of the capacity.\n\nOil dominates the capacity by far, providing almost 70% of the state's."},{"metadata":{"trusted":true},"cell_type":"code","source":"plot = gp_f.groupby('primary_fuel')['capacity_mw'].sum().plot.pie(autopct='%1.0f%%')\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can also take a look at total historical generation, which shows oil was by far the most used source in the past. The current distribution of capacity, with less proportion of oil, shows that the source distribution has been shifting in recent years."},{"metadata":{"trusted":true},"cell_type":"code","source":"plot = gp_f.groupby('primary_fuel')['estimated_generation_gwh'].sum().plot.pie(autopct='%1.0f%%')\nplot.set_title(\"Historical generation by source\")\nplot.set_ylabel(\"\")\npl.pyplot.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can also see that PREPA owns most of the plants in Puerto Rico."},{"metadata":{"trusted":true},"cell_type":"code","source":"gp_f['owner'].value_counts().plot.bar()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"hist = gp_f['capacity_mw'].plot.hist()\nhist.set_xlabel(\"MW\")\npl.pyplot.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Sentinel-5 Precursor NO2 data**"},{"metadata":{},"cell_type":"markdown","source":"This dataset altogether with meteorological data will help us get more insights and calculate emissions (E) and the emissions reduction efficiency (ER) to allow computing emission factors."},{"metadata":{"trusted":true},"cell_type":"code","source":"from datetime import datetime\nimport os\n\nfiles=[]\nfor dirname, _, filenames in os.walk('/kaggle/input/ds4g-environmental-insights-explorer/eie_data/s5p_no2'):\n    for filename in filenames:\n        files.append(os.path.join(dirname, filename))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Here we display band 1 (total vertical NO2 density) of a random file."},{"metadata":{"trusted":true},"cell_type":"code","source":"import rasterio as rio\nfrom matplotlib import pyplot as plt\n\npr_bounds = [[18.6,-67.3,],[17.9,-65.2]]\npr_center = [np.mean(x) for x in zip(*pr_bounds)]\nprint(pr_center)\n\ndef overlay_on_pr(file_name,band_layer):\n    band = rio.open(file_name).read(band_layer)\n    print(band)\n    m = folium.Map(pr_center, zoom_start=8)\n    folium.raster_layers.ImageOverlay(\n        image=band,\n        bounds=pr_bounds,\n        colormap=lambda x: (1, 0, 0, x)\n    ).add_to(m)\n    return m\n\ndef plot_scaled(file_name):\n    vmin, vmax = np.nanpercentile(file_name, (5,95))  # 5-95% stretch\n    img_plt = plt.imshow(file_name, cmap='Purples_r', vmin=vmin, vmax=vmax)\n    plt.show()\n    \nimage_band = rio.open(files[10]).read(1)\nplot_scaled(image_band)\noverlay_on_pr(files[10], 1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from dateutil.parser import parse\n# Process and create DataFrame\nlist_no2 = []\nprint(\"There are\", len(files), \"NO2 Sentinel measurements.\")\nfor file in files:\n    no2 = dict()\n    filename = os.path.splitext(os.path.basename(file))[0]\n    no2['start_ts'] = parse(filename.split('_')[-2])\n    no2['end_ts'] = parse(filename.split('_')[-1]) \n    no2['data'] = rio.open(file).read(1)\n    no2['sum_no2'] = np.sum(no2['data'])\n    no2['avg_no2'] = np.mean(no2['data'])\n    no2['month'] = no2['start_ts'].strftime('%Y-%m')\n    list_no2.append(no2)\nno2_df = pd.DataFrame(list_no2)\nprint(\"Processing files and loading into DF finished.\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.figure(figsize=[8,4])\nchart = sns.lineplot('month', 'avg_no2', data=no2_df)\nplt.setp(chart.get_xticklabels(), rotation=45)\nchart","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"A very simple way to get an average emission factor would be to simply divide the average day pollution in the year by the average capacity produced in one day of that year."},{"metadata":{"trusted":true},"cell_type":"code","source":"no2_df['avg_no2'].mean() / (gp_f['capacity_mw'].sum() / 365)","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}