In [1]:
%matplotlib inline
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
import os
from os import listdir
from os.path import isfile, join
import numpy as np
cwd = os.getcwd()
base_path = join(cwd, '..')
data_path = join(base_path, 'Data storage', 'final state data')
fig_export_path = join(base_path, 'Figures')
idx = pd.IndexSlice
In [2]:
file_date = '2018-03-06'
In [3]:
us_state_abbrev = {
'United States': 'US',
'Alabama': 'AL',
'Alaska': 'AK',
'Arizona': 'AZ',
'Arkansas': 'AR',
'California': 'CA',
'Colorado': 'CO',
'Connecticut': 'CT',
'Delaware': 'DE',
'Florida': 'FL',
'Georgia': 'GA',
'Hawaii': 'HI',
'Idaho': 'ID',
'Illinois': 'IL',
'Indiana': 'IN',
'Iowa': 'IA',
'Kansas': 'KS',
'Kentucky': 'KY',
'Louisiana': 'LA',
'Maine': 'ME',
'Maryland': 'MD',
'Massachusetts': 'MA',
'Michigan': 'MI',
'Minnesota': 'MN',
'Mississippi': 'MS',
'Missouri': 'MO',
'Montana': 'MT',
'Nebraska': 'NE',
'Nevada': 'NV',
'New Hampshire': 'NH',
'New Jersey': 'NJ',
'New Mexico': 'NM',
'New York': 'NY',
'North Carolina': 'NC',
'North Dakota': 'ND',
'Ohio': 'OH',
'Oklahoma': 'OK',
'Oregon': 'OR',
'Pennsylvania': 'PA',
'Rhode Island': 'RI',
'South Carolina': 'SC',
'South Dakota': 'SD',
'Tennessee': 'TN',
'Texas': 'TX',
'Utah': 'UT',
'Vermont': 'VT',
'Virginia': 'VA',
'Washington': 'WA',
'West Virginia': 'WV',
'Wisconsin': 'WI',
'Wyoming': 'WY',
}
In [4]:
index_path = join(data_path, 'Monthly index states {}.csv'.format(file_date))
monthly_index = pd.read_csv(index_path)
In [5]:
annual_index = monthly_index.groupby(['state', 'year']).sum()
annual_index.drop(['month', 'quarter'], axis=1, inplace=True)
annual_index['index (g/kwh)'] = (annual_index['final co2 (kg)']
/ annual_index['generation (mwh)'])
annual_index.head()
Out[5]:
In [6]:
gen_path = join(data_path, 'Monthly generation states {}.csv'.format(file_date))
monthly_gen = pd.read_csv(gen_path)
In [7]:
annual_gen = monthly_gen.groupby(['fuel category', 'state', 'year']).sum()
annual_gen.drop(['month', 'quarter'], axis=1, inplace=True)
annual_gen.head()
Out[7]:
Use data from LBNL to determine what year each state's RPS required that renewable energy be procured.
In [8]:
path = os.path.join(base_path, 'Data storage',
'rps_compliance_data_july_2017.xlsx')
rps = pd.read_excel(path, header=35, usecols='A:V', na_values=['-'])
rps.index = rps.index.droplevel([1, 2])
rps.index.names = ['state', 'Type']
rps_tidy = pd.melt(rps.xs('Total RPS', level='Type').reset_index(),
id_vars='state', var_name='year', value_vars=rps.columns,
value_name='Generation').dropna().sort_values(['state', 'year'])
rps_start = {}
for state in rps_tidy['state'].unique():
first_year = rps_tidy.loc[rps_tidy['state'] == state, 'year'].min()
rps_start[state] = first_year
In [9]:
rps_start
Out[9]:
In [10]:
sns.set()
sns.set_style('white')
This version of the plot is probably better/more flexible
In [124]:
def dumbell_plot(data, years, axis_labels, legend_loc=[], offset_divider=35,
rps_start={}, fig_kwargs={}, figsize=(5,9), legend=True, rps_legend=False,
text_h_align='right', palette='deep', style=None):
'''
This is an example to create a dumbbell chart in Python.
If you would like to provide your data and customize the graph, modify the variables in the section below.
Please be aware that you need matplotlib installed in order for this to work.
'''
prop_cycle = plt.rcParams['axes.prop_cycle']
if style:
plt.style.use(style)
colors = sns.color_palette()
else:
colors = sns.color_palette(palette)
# Styles to be used when plotting the different elements of the graph.
axis_label_style = dict(horizontalalignment=text_h_align,
verticalalignment='center', fontsize=10)
data = data.loc[:, years]
min_data = data.min(axis=1)
max_data = data.max(axis=1)
# Create the figure
fig, ax = plt.subplots(figsize=figsize, **fig_kwargs)
index = range(len(axis_labels))
# Auto-set the state abbr text offset
label_offset = data.max().max() / offset_divider
# Loop N times
for i, (data, year) in enumerate(zip(data.T.values, years)):
color = colors[i]
for value, label, j in zip(data, axis_labels, index):
facecolor = None
if label in rps_start and rps_start[label] <= year:
facecolor = 'w'
ax.scatter(value, j, facecolors=facecolor, zorder=3, color=color,
linewidth=2, s=50)
plt.hlines(y=j, xmin=min_data[j], xmax=max_data[j], zorder=2)
if i == 0:
ax.text(min_data[j] - label_offset, j, label,
**axis_label_style)
plt.yticks(index, ['' for x in axis_labels])
if legend:
for i, year in enumerate(years):
ax.scatter(x=legend_loc[i], y=51, color=colors[i], zorder=3,
s=50, linewidth=2)
plt.text(x=legend_loc[i], y=52, s=str(year), ha='center')
plt.hlines(y=51, xmin=legend_loc[0], zorder=2, xmax=legend_loc[-1])
# Add filled and hollow circles to show RPS status in legend
if rps_legend:
ax.scatter(x=legend_loc[i], y=49, color=colors[i], s=50,
linewidth=2)
plt.text(x=legend_loc[i] * 1.4, y=48.9, s='No RPS', ha='left', va='center')
ax.scatter(x=legend_loc[i], y=47.5, color=colors[i], s=50,
linewidth=2, facecolor='w')
plt.text(x=legend_loc[i] * 1.4, y=47.4, s='RPS', ha='left', va='center')
In [12]:
annual_index.head()
Out[12]:
In [80]:
barbell_index = annual_index.pivot_table(values='index (g/kwh)',
index='state', columns='year')
barbell_index.sort_values(by=[2017], inplace=True)
barbell_index.head()
Out[80]:
In [136]:
len([x for x in rps_start.values() if x <=2008])
Out[136]:
In [135]:
rps_start
Out[135]:
In [133]:
states_index = list(barbell_index.index)
dumbell_plot(barbell_index, [2001, 2008, 2017], states_index,
legend_loc=[500, 300, 100], rps_legend=True,
rps_start=rps_start, palette='colorblind')
plt.vlines(439, -1, 50, colors=['0.5'], zorder=1, #linestyles='dashed',
linewidth=1)
plt.ylim(-1, 53)
# plt.text(x=200, y=40, s='2017\nNational Average\n(439 g CO$_2$ $\mathregular{kWh^{-1}}$)',
# ha='center', va='center', size=11)
plt.text(x=680, y=3, s='2017\nNational Average\n(439 g CO$_2$ $\mathregular{kWh^{-1}}$)',
ha='center', va='center', size=11)
sns.despine(left=True)
plt.xlabel('g CO$_2$ $\mathregular{kWh^{-1}}$')
# label as part of a larger figure
plt.text(-0.02, 0.89, s='a)', ha='right', size=11, transform=ax.transAxes)
path = join(cwd, '..', 'Figures',
'State CO2 intensity {}.pdf'.format(file_date))
plt.savefig(path, bbox_inches='tight')
In [43]:
barbell_index.loc['USA', 2001] = 630
barbell_index.loc['USA', 2008] = 580
barbell_index.loc['USA', 2017] = 439
barbell_index.sort_values(by=[2017], inplace=True)
In [44]:
states_index = list(barbell_index.index)
dumbell_plot(barbell_index, [2001, 2008, 2017], states_index,
legend_loc=[500, 300, 100],
rps_start=rps_start, palette='colorblind')
plt.vlines(439, -1, 50, colors=['0.5'], zorder=1, #linestyles='dashed',
linewidth=1)
plt.ylim(-1, 53)
plt.text(x=200, y=40, s='2017\nNational Average\n(439 g CO$_2$ $\mathregular{kWh^{-1}}$)',
ha='center', size=11)
sns.despine(left=True)
plt.xlabel('g CO$_2$ $\mathregular{kWh^{-1}}$')
path = join(cwd, '..', 'Figures',
'State CO2 intensity {}.pdf'.format('David Hawkins'))
plt.savefig(path, bbox_inches='tight')
In [15]:
barbell_index['change'] = barbell_index[2017] - barbell_index[2001]
barbell_index['% change'] = barbell_index['change'] / barbell_index[2001]
In [16]:
max_change_state = barbell_index.sort_values('change').index[0]
max_change_value = barbell_index.sort_values('change')['change'].values[0]
start = barbell_index.loc[max_change_state, 2001]
end = barbell_index.loc[max_change_state, 2017]
print('The largest absolute reduction is {:0.1f} g/kWh in {}, from {:.1f} to {:.1f}'
.format(max_change_value, max_change_state, start, end))
In [17]:
min_change_state = barbell_index.sort_values('change').index[-1]
min_change_value = barbell_index.sort_values('change')['change'].values[-1]
start = barbell_index.loc[min_change_state, 2001]
end = barbell_index.loc[min_change_state, 2017]
print('The smallest absolute reduction is {:0.1f} g/kWh in {}, from {:.1f} to {:.1f}'
.format(min_change_value, min_change_state, start, end))
In [18]:
max_relchange_state = barbell_index.sort_values('% change').index[0]
max_relchange_value = barbell_index.sort_values('% change')['% change'].values[0]
start = barbell_index.loc[max_relchange_state, 2001]
end = barbell_index.loc[max_relchange_state, 2017]
print('The largest relative reduction is {:0.1%} g/kWh in {}, from {:.1f} to {:.1f}'
.format(max_relchange_value, max_relchange_state, start, end))
In [19]:
min_relchange_state = barbell_index.sort_values('% change').index[-1]
min_relchange_value = barbell_index.sort_values('% change')['% change'].values[-1]
start = barbell_index.loc[min_relchange_state, 2001]
end = barbell_index.loc[min_relchange_state, 2017]
print('The smallest relative reduction is {:0.1%} g/kWh in {}, from {:.1f} to {:.1f}'
.format(min_relchange_value, min_relchange_state, start, end))
In [45]:
path = os.path.join(base_path, 'Data storage', 'Derived data',
'State population.csv')
pop = pd.read_csv(path)
pop.columns = pop.columns.str.lower()
In [46]:
pop.head()
Out[46]:
In [47]:
pop.tail()
Out[47]:
In [48]:
pop['state'] = pop['state'].map(us_state_abbrev)
In [49]:
pop.head()
Out[49]:
In [50]:
annual_index_pop = annual_index.reset_index().merge(pop, on=['state', 'year'])
In [51]:
annual_index_pop.head()
Out[51]:
In [52]:
annual_index_pop['tonne CO2/pop'] = (annual_index_pop['final co2 (kg)']
/ 1000
/ annual_index_pop['population'])
annual_index_pop['MWh/pop'] = (annual_index_pop.loc[:, 'generation (mwh)']
/ annual_index_pop['population'])
In [53]:
annual_index_pop.describe(percentiles=[.1, .25])
Out[53]:
In [101]:
barbell_pop = annual_index_pop.pivot_table(values='tonne CO2/pop',
index='state', columns='year')
barbell_pop.sort_values(by=2017, inplace=True)
barbell_pop.drop(['WY', 'ND', 'WV'], inplace=True)
# states_pop = ['{} '.format(x) for x in barbell_pop.index]
states_pop = list(barbell_pop.index)
# rps_states = list(rps_tidy['state'].unique())
dumbell_plot(barbell_pop, [2001, 2008, 2017], states_pop, offset_divider=30,
rps_start=rps_start, legend=False, palette='colorblind')
plt.ylim(-1, 53)
plt.xlim(None, 25)
sns.despine(left=True)
plt.xlabel('Tonne $\mathregular{CO_2 \ Capita^{-1}}$')
# label as part of a larger figure
plt.text(-0.02, 0.89, s='b)', ha='right', size=11, transform=ax.transAxes)
path = join(fig_export_path,
'State CO2 per capita {}.pdf'.format(file_date))
plt.savefig(path, bbox_inches='tight')
In [77]:
barbell_pop = annual_index_pop.pivot_table(values='tonne CO2/pop',
index='state', columns='year')
barbell_pop.sort_values(by=2017, inplace=True)
barbell_pop = barbell_pop.loc[['WV', 'ND', 'WY']]
states_pop = ['WY', 'ND', 'WV']
# rps_states = list(rps_tidy['state'].unique())
dumbell_plot(barbell_pop, [2001, 2008, 2017], states_pop, offset_divider=30,
rps_start=rps_start, legend=False, palette='colorblind',
figsize=(2.5, 9))
plt.ylim(-1, 53)
# plt.xlim(None, 25)
sns.despine(left=True)
plt.xlabel('Tonne $\mathregular{CO_2 \ Capita^{-1}}$')
path = join(fig_export_path,
'State CO2 per capita insert {}.pdf'.format(file_date))
plt.savefig(path, bbox_inches='tight')
Data is from EIA State Energy Data System. Most recent update was March 9, 2018. Code are described in EIA documentation. Negative values indicate flows out of the state.
In [31]:
path = os.path.join(base_path, 'Data storage', 'State energy flows',
'use_all_phy_update.csv')
flows = pd.read_csv(path)
# ELISP is electricity flows
flows = flows.loc[flows['MSN'] == 'ELISP']
# Drop all years before 2001
flows.drop([str(x) for x in range(1960, 2001)], axis=1, inplace=True)
flows.drop(['Data_Status', 'MSN'], axis=1, inplace=True)
# Electricity flows are given in million kWh (see documentation)
flows.loc[:, '2001':] *= 1000
flows.columns = flows.columns.str.lower()
In [32]:
flows = flows.melt(id_vars='state', var_name='year',
value_name='MWh flow')
flows['year'] = flows['year'].astype(int)
In [33]:
flows.tail()
Out[33]:
In [34]:
annual_index_pop_flows = pd.merge(annual_index_pop, flows, on=['state', 'year'])
annual_index_pop_flows['MWh flow/pop'] = (annual_index_pop_flows['MWh flow']
/ annual_index_pop_flows['population'])
annual_index_pop_flows['share flow'] = (annual_index_pop_flows['MWh/pop']
/ annual_index_pop_flows['MWh flow/pop'])
In [35]:
annual_index_pop_flows.tail()
Out[35]:
In [36]:
barbell_flows = annual_index_pop_flows.pivot_table(values='MWh flow/pop',
index='state', columns='year')
barbell_flows.head()
Out[36]:
In [37]:
barbell_flows = annual_index_pop_flows.pivot_table(values='MWh flow/pop',
index='state', columns='year')
barbell_flows.sort_values(by=2015, ascending=False, inplace=True)
# barbell_flows.drop('VT', inplace=True)
states_flows = list(barbell_flows.index)
dumbell_plot(barbell_flows, [2001, 2008, 2015], states_flows, offset_divider=8,
rps_start=rps_start, legend=True, palette='colorblind',
legend_loc=[-20, -40, -60])
plt.ylim(-1, 53)
sns.despine(left=True)
plt.xlabel('MWh flow/Capita')
plt.vlines(0, -1, 53, colors=['0.5'], linewidth=1)
# plt.xticks()
# ax = plt.gca()
# ax.set_xscale("log", nonposx='clip')
path = join(fig_export_path, 'SI', 'State MWh flow per capita.pdf')
# plt.savefig(path, bbox_inches='tight')