{"cells":[{"metadata":{},"cell_type":"markdown","source":"# SEDAC Introduction - How to Work with this Dataset\n---\nThis dataset contains a number of useful data sources in tabular and geospatial format. Some of the data requires processing and potentially unfamiliar functions for access. This notebook contains guides to get you started for each of the datasets contained here. In order these are:\n\n1. 1x1 Gridded Emissions\n2. Trends in Global Freshwater Availability\n3. Development Threat Index\n4. Environmental Performance Index\n5. Climate Effects on Food Supply\n6. Human Modification of Terrestrial Systems (large TIFF file)\n7. Global Annual PM2.5 Grids (large TIFF file)\n8. Natural Resource Protection & Child Health Indicators\n\nThe large TIFF files can cause RAM crashes if you repeatedly load them in, so make sure to only load them once in your code before processing."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import os\nimport math\nimport numpy as np \nimport pandas as pd \nimport matplotlib.pyplot as plt\nimport rasterio #for working with TIFF files\nfrom scipy import interpolate #in case TIFF grid is missing local values, we can impute them","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 1x1 Gridded Emissions\n---\n\nThis dataset consists of ASCII files that create a grid of the world map, where the values in each location correspond to emission levels in that area. This is also accompanied by two CSV files that give some overall summary values for different regions of the world, and their total predicted emissions for various pollutants over sequential decades for different SRES scenarios. For more information on these projected scenarios, please consult the [IPCC guide (PDF)](https://www.ipcc.ch/site/assets/uploads/2018/03/emissions_scenarios-1.pdf).\n\nThe CSV files must first be processed as they are in a non-standard format. The desired values can then be selected using ordinary indexing:"},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"def process_so2(df):\n    cols = ['region', 'scenario', '1990', '2000', '2010', '2020', '2030', '2040', '2050', '2060', '2070', '2080','2090', '2100']\n    df.columns = cols\n    df[['region', 'scenario']] = df[['region', 'scenario']].replace(math.nan, 'nan')\n    for i in range(len(df)):\n        scenario = df.iloc[i, :].scenario\n        if scenario!='nan':\n            current_scenario = scenario\n        else:\n            df.loc[i, 'scenario'] = current_scenario\n    return df.loc[df.region!='nan', :]\n\n\ndef process_rg(df):\n    years = np.arange(1990, 2110, 10).astype(str).tolist()\n    df[['Scenario']] = df[['Scenario']].replace(math.nan, 'nan')\n    row_list = []\n    for i in range(len(df)):\n        scenario = df.iloc[i, :]['Scenario']\n        gas_region = df.iloc[i, :]['Gas/Region']\n        readings = df.iloc[i, :][years].values\n        if scenario!='nan':\n            current_scenario = scenario\n            current_gas = gas_region\n            row = pd.DataFrame(index=np.arange(1), columns=['scenario', 'region', 'gas'] + years)\n            row['scenario'] = current_scenario\n            row['gas'] = current_gas\n            row['region'] = np.nan\n        else:\n            row = pd.DataFrame(index=np.arange(1), columns=['scenario', 'region', 'gas'] + years)\n            row['scenario'] = current_scenario\n            row['gas'] = current_gas\n            row['region'] = gas_region\n            row[years] = readings\n            \n        row_list.append(row)\n                        \n    out_df = pd.concat(row_list) \n    out_df[years] = out_df[years].astype('float32')\n    return out_df.loc[~out_df.region.isna(), :]\n\nso2_df = pd.read_csv('/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/1x1 Gridded Emissions/so2-values.csv')\nreactive_gas_df = pd.read_csv('/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/1x1 Gridded Emissions/reactive-gases.csv')\n\nso2 = process_so2(so2_df)\nreactive_gas_df = process_rg(reactive_gas_df)\nso2.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reactive_gas_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"markdown","source":"### Working with ASCII data\n---\nFirst, for convenience we can make a dataframe from the available files:"},{"metadata":{"trusted":true},"cell_type":"code","source":"scenarios = ['A1AIM', 'A2ASF', 'B2MESSAGE', 'A1GMINICAM', 'B1IMAGE', 'A1TMESSAGE']\ndict_list = []\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        for scenario in scenarios:\n            if filename[:len(scenario)]==scenario:\n                temp_dict = {}\n                temp_dict['scenario'] = scenario\n                temp_dict['year'] = filename[len(scenario):len(scenario)+4]\n                temp_dict['gas'] = filename[len(scenario)+4:].split('.')[0]\n                temp_dict['path'] = os.path.join(dirname, filename)\n                dict_list.append(temp_dict)\n                \ngrid_df = pd.DataFrame(dict_list)\ngrid_df                ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now we can get the data for a particular year, scenario and pollutant:"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_emissions(path):\n    out_df = pd.DataFrame(np.loadtxt(grid_df.iloc[0, :].path, delimiter=',', skiprows=12), columns=['long', 'lat', 'value'])\n    out_df['long'] = out_df['long'].astype(int)\n    out_df['lat'] = out_df['lat'].astype(int)\n    return out_df\n    \ndef get_path(grid_df, year, scenario, gas):\n    return grid_df.loc[(grid_df.year==year) & (grid_df.scenario==scenario) & (grid_df.gas==gas), :].path\n    \ntest_path = get_path(grid_df, 2080, 'A2ASF', 'NVMOC')\na2asf_2080_nvmoc_df = get_emissions(test_path)\na2asf_2080_nvmoc_df.head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"These can now be matched with any relevant coordinates for cities and corporate facilities.\n\n\nWe can also visualize the data to ensure it looks correct:"},{"metadata":{"trusted":true},"cell_type":"code","source":"def show_grid(df):\n    df['long'] = df['long'] - df['long'].min()\n    df['lat'] = df['lat'] - df['lat'].min()\n    x_range = np.arange(0, df.long.max()+1)\n    y_range = np.arange(0, df.lat.max()+1)\n    grid = np.zeros([len(y_range), len(x_range)])\n    print(grid.shape)\n    for x in np.unique(df['long'].values):\n        temp_df = df.loc[df.long==x, :]\n        lats = temp_df['lat']\n        values = temp_df['value']\n        grid[lats, x] = values\n    \n    grid = np.flipud(grid)\n    plt.figure(figsize=(15, 6))\n    plt.imshow(grid)\n    plt.show();\n    \nshow_grid(a2asf_2080_nvmoc_df)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"markdown","source":"This particular pollutant has high predicted concentrations across large regions of Asia, with spikes in North America, Europe and Sub-Saharan Africa.\n\n# Trends in Global Freshwater Availability\n---\nA GeoTIFF is a TIFF image file that has been augmented with some additional tags which denote properties such as coordinate system and geographical extent. The rasterio package can be used to load in these files. Since the actual values are stored in a NumPy array, accessing and manipulating them is straightforward:"},{"metadata":{"trusted":true},"cell_type":"code","source":"fp = '/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Trends in Global Freshwater Availability/freshwater_availability.tif'\nimg = rasterio.open(fp)\nimg_array = img.read(1)\nprint(f'Maximum change: {img_array.max()}')\nprint(f'Minimum change: {img_array.min()}')\nprint(f'Raster shape: {img_array.shape}')\nfig = plt.figure(figsize=(15, 6))\nax = fig.add_subplot()\nim = ax.imshow(img_array)\nfig.colorbar(im)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Visually, this figure is somewhat dominated by the large projected decline in freshwater reserves at the polar regions. The researchers recommend using a mask to remove all ocean and sea coverage, as the instruments used for these measurements detect alterations in the local gravity caused by water flows, which is naturally inaccurate over moving bodies of water.\n\nHowever, this data can still yield useful insights into the rate of freshwater decline in a particular region. Looking at the shape, it is a simple array of dimensions `(360, 720)` which allow values to be mapped to coordinates for cities and other locations. For this function, longitudes West must be made negative. Say, for Chicago at 41.8781° N, 87.6298° W:"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_value(raster, array, lon, lat):\n    # get numpy index of these coordinates\n    py, px = raster.index(lon, lat)\n    return array[py, px]\n\nget_value(img, img_array, -87.6298, 41.8781)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Development Threat Index\n\n---\nThis data is also stored in GeoTIFF format. For this demonstration, we will only consider the cumulative DTI score. More specific ones relating to fossil fuels, mining, agricultural expansion and more can also be found in this directory. They are accompanied by PDF files for the methodology and documentation for each considered variable.\n\nThe values for ocean/sea coverage can be replaced with `np.nan` for better visibility:"},{"metadata":{"trusted":true},"cell_type":"code","source":"fp = '/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Development Threat Index/geotiff/lulc-development-threat-index_geogrpahic.tif'\nimg = rasterio.open(fp)\nimg_array = img.read(1).astype('float32')\nimg_array[img_array==128] = np.nan\nprint(f'Raster shape: {img_array.shape}')\nfig = plt.figure(figsize=(15, 6))\nax = fig.add_subplot()\nim = ax.imshow(img_array)\nfig.colorbar(im)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"There are still some artifacts present in the Pacific ocean, probably a consequence of the transformation from the original Molleweide projection. Despite being an unusual shape, it is simple to map coordinates to their relevant Development Threat Index score. Again, for Chicago:"},{"metadata":{"trusted":true},"cell_type":"code","source":"get_value(img, img_array, -87.6298, 41.8781)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Environmental Performance Index\n---\n\nThe EPI data is stored in .xlsx format, combining multiple spreadsheets into a single file. This can best be accessed using `pd.read_excel`, making sure to use the parameter `sheet_name=None`. This will construct a dictionary of dataframes that can be accessed via the dictionary key:"},{"metadata":{"trusted":true},"cell_type":"code","source":"epi_d = pd.read_excel('/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Environmental Performance Index/2018-epi.xlsx', sheet_name=None)\nepi_d.keys()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now if we want to examine the `'2018EPI_ScoresCurrent'` table:"},{"metadata":{"trusted":true},"cell_type":"code","source":"epi_d['2018EPI_ScoresCurrent']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Details on the different fields can be found in the accompanying PDF, or examined directly in the time-series data sheet. The latter also contains a wealth of data:"},{"metadata":{"trusted":true},"cell_type":"code","source":"epi_ts_d = pd.read_excel('/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Environmental Performance Index/2018-epi-raw-data-time-series.xlsx', sheet_name=None)\nprint(epi_ts_d.keys())\nepi_ts_d['DataDictionary']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"So if we want to examine the Wastewater Treatment/Water Resources score for the US:"},{"metadata":{"trusted":true},"cell_type":"code","source":"epi_df = epi_d['2018EPI_ScoresCurrent']\nepi_df.loc[epi_df.iso=='USA']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can see that the Wastewater Treatment score is 92.57. Its total Environmental  Performance Index score meanwhile is 71.19. There are many features that can be extracted from this dataset, and as it is dated from 2018 it can be considered an up-to-date source.\n\n# Climate Effects on Food Supply\n---\nThis dataset examines the projected changes in the supply of essential crops worlwide and by country. This data is a bit older than the others and is based on the SRES scenarios created for IPCC projections. You can however create useful features by examining the difference between these predictions and known crop yields in recent years - [Earthstat](http://www.earthstat.org/) would be a good place to start.\n\nIt is also in .xlsx format - please note the spaces in front of the table names `' Data Dictionary'` and `' Data sheet'`:"},{"metadata":{"trusted":true},"cell_type":"code","source":"fs_d = pd.read_excel('/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Climate Effects on Food Supply/crops-yield-changes-hadcm3-sres.xls', sheet_name=None)\nprint(fs_d.keys())\nfs_d[' Data Dictionary']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"fs_d[' Data sheet']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"So, if we want to predict the weighted average yield change for grain crops in the US by 2020:"},{"metadata":{"trusted":true},"cell_type":"code","source":"fs_df = fs_d[' Data sheet']\nfs_df.loc[fs_df['BLS_2_Countries_(SRES)_ABBREVNAME']=='United States']['B2B2020']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can see that it was expected to be around 18% higher than the known values for 2000.\n\n# Human Modification of Terrestrial Systems\n---\nThis dataset scores areas of the Earth's land surface between 0 and 1, a combined metric that includes the effects of 13 identified 'environmental stressors'. It is a useful score for determining the local effect on the environment caused by city or corporate activity. The dataset presents some challenges since not everywhere is scored, so we have to interpolate some results."},{"metadata":{"trusted":true},"cell_type":"code","source":"fp = '/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Human Modification of Terrestrial Sytems/lulc-human-modification-terrestrial-systems_geographic.tif'\nimg = rasterio.open(fp)\nimg_array = img.read(1).astype('float32')\nimg_array[img_array<0] = np.nan #oceanic areas have value -3.4028235e+38\nprint(img_array.shape)\nfig = plt.figure(figsize=(15, 6))\nax = fig.add_subplot()\nim = ax.imshow(img_array)\nfig.colorbar(im)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_hmts(raster, array, long, lat):\n    \n    #if value is present at these coords, return value\n    out = get_value(img, img_array, long, lat)\n    if not np.isnan(out):\n        return \n    \n    #otherwise perform 2d linear interpolation on local area and impute the value\n    py, px = raster.index(lon, lat) \n    array = array[py-100:py+100, px-100:px+100]\n    x = np.arange(0, array.shape[1])\n    y = np.arange(0, array.shape[0])\n    array = np.ma.masked_invalid(array)\n    xx, yy = np.meshgrid(x, y)\n    #get only the valid values\n    x1 = xx[~array.mask]\n    y1 = yy[~array.mask]\n    newarr = array[~array.mask]\n\n    interpolated_arr = interpolate.griddata((x1, y1),\n                                            newarr.ravel(),\n                                            (xx, yy),\n                                            method='linear')            \n    #return central value from the imputed local area\n    return interpolated_arr[100, 100]\n\nget_hmts(img, img_array, -87.6298, 41.8781)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Global Annual PM2.5 Grids\n---\nPM2.5 refers to the presence of fine, particulate matter either natural or man-made which can cause health problems in high concentrations/long exposures. Specifically, PM2.5 is the mass per cubic metre of air of particles with a diameter less than 2.5 micrometres (µm). This dataset supplies a grid for these readings, indicating a valuable property of general air quality."},{"metadata":{"trusted":true},"cell_type":"code","source":"fp = '/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Global Annual PM2.5 Grids/sdei-global-annual-gwr-pm2-5-modis-misr-seawifs-aod-2016-geotiff/gwr_pm25_2016.tif'\nimg = rasterio.open(fp)\nimg_array = img.read(1).astype('float32')\nimg_array[img_array<-1e6] = np.nan #oceanic areas have value -3.4028235e+38\nprint(img_array.shape)\nfig = plt.figure(figsize=(15, 6))\nax = fig.add_subplot()\nim = ax.imshow(img_array)\nfig.colorbar(im)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Similarly to the Human Modification of Terrestrial Systems dataset, there are some missing values over land areas so we can interpolate values if they are missing from the original raster:"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_pm25(raster, array, long, lat):\n    \n    #if value is present at these coords, return value\n    out = get_value(img, img_array, long, lat)\n    if not np.isnan(out):\n        return \n    \n    #otherwise perform 2d linear interpolation on local area and impute the value\n    py, px = raster.index(lon, lat) \n    array = array[py-100:py+100, px-100:px+100]\n    x = np.arange(0, array.shape[1])\n    y = np.arange(0, array.shape[0])\n    array = np.ma.masked_invalid(array)\n    xx, yy = np.meshgrid(x, y)\n    #get only the valid values\n    x1 = xx[~array.mask]\n    y1 = yy[~array.mask]\n    newarr = array[~array.mask]\n\n    interpolated_arr = interpolate.griddata((x1, y1),\n                                            newarr.ravel(),\n                                            (xx, yy),\n                                            method='linear')            \n    #return central value from the imputed local area\n    return interpolated_arr[100, 100]\n\nget_pm25(img, img_array, -87.6298, 41.8781)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Natural Resource Protection & Child Health Indicators\n---\nThe final dataset in this collection provides scoring metrics for natural resource protection and child health indicators. The dataset was uploaded this year with the latest complete values being for 2019, making it an up-to-date resource. \n\nIt is stored in .xlsx format and can be accessed like the others:"},{"metadata":{"trusted":true},"cell_type":"code","source":"nrp_chi_d = pd.read_excel('/kaggle/input/socioeconomic-data-and-applications-center/sedac-data/Natural Resource Protection and Child Health Indicators/nrpi-chi-2019-xlsx.xlsx', sheet_name=None)\nprint(nrp_chi_d.keys())\nnrp_chi_d['data dictionary']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The 2019 scores for NRPI and CHI can be trivially accessed by country:"},{"metadata":{"trusted":true},"cell_type":"code","source":"nrp_chi_d['NRPI 2019'].head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"nrp_chi_d['CHI 2019'].head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"In these columns, the 'v2019' refers to the dataset release, while the last two digits are for the year in question. So the latest values for NRPI and CHI for the United States can be extracted as follows:"},{"metadata":{"trusted":true},"cell_type":"code","source":"nrp_df = nrp_chi_d['NRPI 2019']\nchi_df = nrp_chi_d['CHI 2019']\nprint(f\"United states latest Natural Resource Protection Indicator: {nrp_df.loc[nrp_df.ISO3=='USA', :].NRPI_v2019_19.values[0]}\")\nprint(f\"United states latest Child Health Indicator: {chi_df.loc[chi_df.ISO3=='USA', :].CHI_v2019_18.values[0]}\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"A full guide to the columns can be found in the accompanying PDF, and in the `'data dictionary'` dataframe.\n\n# Conclusion\n---\n\nThis collection of SEDAC resources contains a wealth of information about socioeconomic and environmental data worldwide, including current metrics and projected future values. The geospatial data can be used to enrich both visualizations and tabular data, when it is matched to coordinates using the functions provided in this notebook. Each dataset contains an explanatory PDF in case you wish to read further about the methods used for the data collection and generation. \n\nI hope this dataset is of use to you in your work. If you have any comments or questions, please leave them in the section below."}],"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}