notebooks/Data_Prep.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "from itertools import combinations\n",
    "import datetime as dt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def drop_outliers(df, columns, frac = 0.05):\n",
    "    \n",
    "    df = df.copy()\n",
    "    \n",
    "    # Mark outliers for each column\n",
    "    for column in columns:\n",
    "        cut_low, cut_high = df[column].quantile([frac/2, (1-frac/2)])\n",
    "        df[column + '_keep'] = ((df[column] < cut_high) & (df[column] > cut_low)).astype(int)\n",
    "        \n",
    "    # Subset to observations that aren't outliers for any column\n",
    "    query_str = ' & '.join([x + '_keep == 1' for x in columns])\n",
    "    \n",
    "    return df.query(query_str).drop([x + '_keep' for x in columns], axis = 1).copy()\n",
    "\n",
    "def df_drop_na(df, cols):\n",
    "    \n",
    "    return df[~df[cols].isna().any(axis=1)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Import Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "price_df       = pd.read_csv('../data/electricity/elc_retail_price_monthly.csv', na_values = 'NM')\n",
    "sales_df       = pd.read_csv('../data/electricity/elc_retail_sales_monthly.csv', na_values = 'NM')\n",
    "fossil_cost_df = pd.read_csv('../data/instruments/fossil_fuel_cost_monthly.csv', na_values = ['W', '--', 'NM'])\n",
    "pcepi_df       = pd.read_csv('../data/other/us_pcepi.csv')\n",
    "ddd_df         = pd.read_csv('../data/controls/state_monthly_hdd_cdd.csv')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clean and Merge Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fix state names\n",
    "sales_df['State']       = sales_df['State'].apply(lambda x: x.strip().lower()) \n",
    "price_df['State']       = price_df['State'].apply(lambda x: x.strip().lower()) \n",
    "fossil_cost_df['state'] = fossil_cost_df['state'].apply(lambda x: x.strip().lower()) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fix dates\n",
    "pcepi_df['month'] = pcepi_df['DATE'].apply(lambda x: x.split('/')[0])\n",
    "pcepi_df['year']  = pcepi_df['DATE'].apply(lambda x: x[-4:])\n",
    "ddd_df['month']   = ddd_df['month'].astype(str)\n",
    "ddd_df['year']    = ddd_df['year'].astype(str)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Reshape dataframes\n",
    "sales_melt_df  = pd.melt(sales_df, id_vars = ['State'], value_vars = sales_df.columns[1:],\n",
    "                         var_name = 'date', value_name = 'load')\n",
    "\n",
    "\n",
    "price_melt_df  = pd.melt(price_df, id_vars = ['State'], value_vars = price_df.columns[1:],\n",
    "                         var_name = 'date', value_name = 'price')\n",
    "\n",
    "\n",
    "fossil_cost_df         = pd.melt(fossil_cost_df, id_vars = ['state', 'fuel'], value_vars = fossil_cost_df.columns[4:],\n",
    "                         var_name = 'date', value_name = 'cost')\n",
    "fossil_cost_df         = fossil_cost_df.reset_index(drop = True)\n",
    "fossil_cost_df['cost'] = fossil_cost_df['cost'].fillna(-1234)\n",
    "fossil_cost_df         = fossil_cost_df.groupby(['state', 'date', 'fuel'])['cost'].sum().unstack('fuel').reset_index()\n",
    "fossil_cost_df.columns = (list(fossil_cost_df.columns[:2]) + \n",
    "                          [x.strip().replace(' ', '_') + '_cost' for x in fossil_cost_df.columns[2:]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Merge electricity and instrument data\n",
    "data_df = price_melt_df.merge(sales_melt_df, on = ['State', 'date']).rename(columns = {'State': 'state'}).merge(\n",
    "            fossil_cost_df, on = ['state', 'date'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Month abbreviations\n",
    "month_abrv_dict = dict(zip(['Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec'],\n",
    "                          range(1,13)))\n",
    "\n",
    "# Clean up columns\n",
    "data_df['state'] = data_df['state'].apply(lambda x: x.strip().lower()) \n",
    "data_df['month'] = data_df['date'].apply(lambda x: str(month_abrv_dict.get(x.split('-')[0])))\n",
    "data_df['year']  = data_df['date'].apply(lambda x: str(2000+int(x.split('-')[1])))\n",
    "data_df['price'] = pd.to_numeric(data_df['price'])/100\n",
    "data_df['load']  = pd.to_numeric(data_df['load'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Merge with pcepi and ddd data\n",
    "data_df = data_df.merge(pcepi_df, on =['month', 'year']).drop('DATE', axis = 1).merge(\n",
    "                        ddd_df, on = ['state', 'month', 'year'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fix missing entries in fuel cost\n",
    "fuel_cost_cols = [x for x in data_df.columns if 'cost' in x]\n",
    "for fuel_cost_col in fuel_cost_cols:\n",
    "    data_df[fuel_cost_col] = pd.to_numeric(fossil_cost_df[fuel_cost_col]).apply(\n",
    "        lambda x: np.nan if (x == -1234 or x == 0) else x)\n",
    "    \n",
    "# Adjust for inflation\n",
    "for price_col in (fuel_cost_cols + ['price']):\n",
    "    data_df[price_col] = data_df[price_col]*100/data_df['PCEPI']\n",
    "    \n",
    "# Fix date columns\n",
    "data_df['month'] = data_df['month'].astype(int) \n",
    "data_df['year']  = data_df['year'].astype(int)  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Construct Regression Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Preliminary Data Set Observations:   818\n"
     ]
    }
   ],
   "source": [
    "# Subset dataframe for relevant variables\n",
    "em_sample_df = data_df[['state', 'price', 'load', 'month', 'year', 'coal_cost', 'CDD', 'HDD']].copy()\n",
    "\n",
    "# Drop missing and outliers (CDD and HDD data already cleaned)\n",
    "em_sample_df[~em_sample_df.isna().any(axis = 1)]\n",
    "em_sample_df = drop_outliers(em_sample_df, columns = ['price', 'load', 'coal_cost'], frac = 0.01).reset_index(drop = True)\n",
    "\n",
    "# Add date column\n",
    "em_sample_df['date'] = em_sample_df.apply(lambda x: '{0}/{1}'.format(x.month, x.year), axis = 1)\n",
    "\n",
    "print('Preliminary Data Set Observations:  ', len(em_sample_df))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Regression Data Set Observations:   6817\n"
     ]
    }
   ],
   "source": [
    "# Get all unique combinations of (t,s,i,j) => (month, month, state, state)\n",
    "sample_array = np.array([np.array(comb) for comb in combinations(em_sample_df.index, 2)])\n",
    "\n",
    "# Split into two dataframes for (t,i) and (s,j)\n",
    "em_sample_df_1         = em_sample_df.loc[sample_array[:,0]].reset_index(drop = True)\n",
    "em_sample_df_2         = em_sample_df.loc[sample_array[:,1]].reset_index(drop = True)\n",
    "em_sample_df_1.columns = [x + '_1' for x in em_sample_df_1.columns] \n",
    "em_sample_df_2.columns = [x + '_2' for x in em_sample_df_2.columns] \n",
    "\n",
    "# Merge and filter out t == s and i != j\n",
    "em_concat_df = pd.concat([em_sample_df_1, em_sample_df_2], axis = 1)\n",
    "em_concat_df = em_concat_df.query('state_1 == state_2 & date_1 != date_2').drop_duplicates().reset_index(drop = True)\n",
    "\n",
    "# Add reg data columns\n",
    "em_concat_df['ln_load_rel']    = np.log(np.divide(em_concat_df['load_1'],      em_concat_df['load_2']))\n",
    "em_concat_df['ln_price_rel']   = np.log(np.divide(em_concat_df['price_1'],     em_concat_df['price_2']))\n",
    "em_concat_df['ln_coal_rel']    = np.log(np.divide(em_concat_df['coal_cost_1'], em_concat_df['coal_cost_2']))\n",
    "\n",
    "# Regressions should be done with pairs of data for the same state\n",
    "reg_data_df = em_concat_df.copy()\n",
    "reg_data_df['time_diff'] = (reg_data_df\n",
    "                            .apply(lambda x: \n",
    "                                   dt.datetime(x.year_1, x.month_1, 1) - dt.datetime(x.year_2, x.month_2, 1), \n",
    "                                   axis = 1)\n",
    "                            .apply(lambda x: np.round(x.days/30.4, 0)))\n",
    "\n",
    "print('Regression Data Set Observations:  ', len(reg_data_df))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "reg_data_df.to_csv('../data/processed/regression_data.csv', index = False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data Description"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>price</th>\n",
       "      <th>load</th>\n",
       "      <th>coal_cost</th>\n",
       "      <th>CDD</th>\n",
       "      <th>HDD</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>848.000000</td>\n",
       "      <td>848.00000</td>\n",
       "      <td>848.000000</td>\n",
       "      <td>848.000000</td>\n",
       "      <td>848.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>0.118393</td>\n",
       "      <td>2454.96934</td>\n",
       "      <td>46.004220</td>\n",
       "      <td>64.662736</td>\n",
       "      <td>443.199292</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>0.027906</td>\n",
       "      <td>2380.24399</td>\n",
       "      <td>17.683376</td>\n",
       "      <td>134.965191</td>\n",
       "      <td>412.267006</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>0.076601</td>\n",
       "      <td>144.00000</td>\n",
       "      <td>24.038135</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>0.098699</td>\n",
       "      <td>790.50000</td>\n",
       "      <td>33.429777</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>48.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>0.110991</td>\n",
       "      <td>1869.50000</td>\n",
       "      <td>41.396633</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>367.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>0.128693</td>\n",
       "      <td>3315.00000</td>\n",
       "      <td>51.831937</td>\n",
       "      <td>47.250000</td>\n",
       "      <td>754.500000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>0.208995</td>\n",
       "      <td>18621.00000</td>\n",
       "      <td>107.489340</td>\n",
       "      <td>761.000000</td>\n",
       "      <td>1794.000000</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "            price         load   coal_cost         CDD          HDD\n",
       "count  848.000000    848.00000  848.000000  848.000000   848.000000\n",
       "mean     0.118393   2454.96934   46.004220   64.662736   443.199292\n",
       "std      0.027906   2380.24399   17.683376  134.965191   412.267006\n",
       "min      0.076601    144.00000   24.038135    0.000000     0.000000\n",
       "25%      0.098699    790.50000   33.429777    0.000000    48.000000\n",
       "50%      0.110991   1869.50000   41.396633    0.000000   367.000000\n",
       "75%      0.128693   3315.00000   51.831937   47.250000   754.500000\n",
       "max      0.208995  18621.00000  107.489340  761.000000  1794.000000"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data_df.query('coal_cost == coal_cost')[['price', 'load', 'coal_cost', 'CDD', 'HDD']].describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\\begin{tabular}{lrrrrrrr}\n",
      "\\toprule\n",
      "{} &       mean &        std &       min &       25\\% &        50\\% &        75\\% &         max \\\\\n",
      "\\midrule\n",
      "price     &     0.1184 &     0.0279 &    0.0766 &    0.0987 &     0.1110 &     0.1287 &      0.2090 \\\\\n",
      "load      &  2454.9693 &  2380.2440 &  144.0000 &  790.5000 &  1869.5000 &  3315.0000 &  18621.0000 \\\\\n",
      "coal\\_cost &    46.0042 &    17.6834 &   24.0381 &   33.4298 &    41.3966 &    51.8319 &    107.4893 \\\\\n",
      "CDD       &    64.6627 &   134.9652 &    0.0000 &    0.0000 &     0.0000 &    47.2500 &    761.0000 \\\\\n",
      "HDD       &   443.1993 &   412.2670 &    0.0000 &   48.0000 &   367.0000 &   754.5000 &   1794.0000 \\\\\n",
      "\\bottomrule\n",
      "\\end{tabular}\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(data_df.query('coal_cost == coal_cost')[['price', 'load', 'coal_cost', 'CDD', 'HDD']]\n",
    "      .describe()\n",
    "      .drop('count', axis = 0)\n",
    "      .T\n",
    "      .round(decimals = 4)\n",
    "      .to_latex())"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}