Project 2: Visualizing Power Plant Emissions vs. Power Plant Generation

Project 2: Visualizing Power Plant Emissions vs. Power Plant Generation#

In this project, I used pandas and plotly to explore the relationship between power plant electricity generation and greenhouse gas emissions according to publicly available datasets from the U.S. Environmental Protection Agency (EPA) and the U.S. Energy Information Administration (EIA). Since both agencies group power plants using a common ID system, I was able to create a merged dataset linking each power plant’s emissions with its net electricity generation.

Then, by using a library to re-categorize each of the EIA’s hyper-specific fuel codes according to commonly understood umbrella terms, like coal, oil, gas, biomass, and renewables, I showed in a scatter plot how a plant’s emissions compared to its electricity production across different fuel types.

On this page, I will walk users step-by-step through the process so that they can easily replicate it if they wish!

First, I imported the packages I planned to use for the data cleaning, merge, and visualization, plotly and pandas.

Pandas is a Python library used for data cleaning, manipulation, and analysis using its core data structure, the DataFrame.

Plotly is an open-source, interactive graphing library used to create publication-quality visualizations like charts and maps from data.

import pandas as pd 
import plotly.express as px

Next, I downloaded from EPA GHGRP’s website its most recent FLIGHT dataset, which provides facility level information on greenhouse gas emissions for polluting facilities across the U.S. I chose the 2023 dataset because it was the most recent available. Because it downloaded as an Excel file, I had to open and resave it as a CSV for ease of use with pandas.

Then, I asked Python and Pandas to skip the first three rows, which were headers with overarching information about the dataset, and so not useful to the task at hand.

As a result, ghg is now a Pandas DataFrame containing all the emissions data from the CSV, with proper column headers starting at row 4 of the original file.

ghg = pd.read_csv("/Users/ariakovalovich/Documents/SIPA MIA Coursework/aak2257" 
".github.io/Project 2 Datasets/ghgp_data_2023.csv", skiprows=3)

print("GHGRP columns:", ghg.columns.tolist()[:12])

GHGRP columns: ['Facility Id', 'FRS Id', 'Facility Name', 'City', 'State', 'Zip Code', 'Address', 'County', 'Latitude', 'Longitude', 'Primary NAICS Code', 'Industry Type (subparts)']

/var/folders/pl/4_5qmb416_v7ywlr51mndq9h0000gn/T/ipykernel_41382/3037226414.py:1: DtypeWarning:

Columns (2,3,4,6,7,11,12,63,64,65) have mixed types. Specify dtype option on import or set low_memory=False.

Next, because I only need a few parameters for my analysis, after loading the full GHGRP dataset, I created a smaller dataframe with only the following columns:

Facility Id: the unique identifier for each reporting facility
Facility Name: the facility’s name
State: the U.S. state where the facility is located
Total reported direct emissions: the total greenhouse gases directly emitted by the facility

I extracted these columns into a smaller dataframe called ghg_small for simpler analysis and easier merging with other datasets. To make the column names consistent, I renamee “Facility Id” to “Facility_Id”, replacing the space with an underscore.

ghg_small = ghg[
    [
        "Facility Id",
        "Facility Name",
        "State",
        "Total reported direct emissions",
    ]
]
ghg_small = ghg_small.rename(columns={"Facility Id": "Facility_Id"})

Through some Googling, I learned that in order to successfully merge EPA and EIA datasets, I needed to add another column to the GHGRP dataset that includes another unique identifier code for each power plant, known as the ORIS Code in EPA datasets and Plant ID in EIA datasets. Luckily, EPA provides ORIS codes in a separate crosswalk dataset, which I downloaded and converted to CSV.

In the below code, I asked Pandas and Python to read the csv and print the column names into a Python list so that i could find the exact name of the relevant Oris Code column.

cross = pd.read_csv(
    "/Users/ariakovalovich/Documents/SIPA MIA Coursework/aak2257.github.io"
    "/Project 2 Datasets/ghgrp_oris_power_plant_crosswalk_12_13_21.csv"
)

print("Crosswalk columns:", cross.columns.tolist())

Crosswalk columns: ['GHGRP Facility ID', 'FACILITY NAME', 'GHGRP - City', 'GHGRP - State', 'ORIS CODE', 'ORIS CODE 2', 'ORIS CODE 3', 'ORIS CODE 4', 'ORIS CODE 5', 'GHGRP - Power Plant Sector', 'Reported Subparts RY10', 'Reported Subparts RY11', 'Reported Subparts RY12', 'Reported Subparts RY13', 'Reported Subparts RY14', 'Reported Subparts RY15', 'Reported Subparts RY16', 'Reported Subparts RY17', 'Reported Subparts RY18', 'Reported Subparts RY19', 'Reported Subparts RY20']

The full crosswalk has many columns, but for merging I only needed “GHGRP Facility ID” to match the GHGRP emissions data and “ORIS CODE” to match the EIA power generation data. The below code creates a smaller dataframe called cross_small that is easier to work with. Then, for consistency, I renamed the GHGRP Facility ID column to Facility_Id to match the ghg_small dataset. I then confirmed that the prior operations worked by confirming the column names and looking at the head of the dataset.

cross_small = cross[["GHGRP Facility ID", "ORIS CODE"]]

cross_small = cross_small.rename(columns={"GHGRP Facility ID": "Facility_Id"})
cross_small.columns
cross_small.head()

	Facility_Id	ORIS CODE
0	1000001	54537
1	1000019	57937
2	1000031	57403
3	1000033	50300
4	1000034	10861

Before I merged the datasets, I checked to make sure that the facility identifier columns in both my GHGRP emissions data (ghg_small) and the crosswalk (cross_small) were in compatible formats. Both datasets use a column called Facility_Id, and I confirmed their data types. This step helped ensure that the merge would work correctly without type-related issues.

ghg_small.Facility_Id.dtype
cross_small.Facility_Id.dtype

dtype('int64')

I then performed an inner merge of the two datasets. By merging on “Facility_Id”, I was able to combine the emissions information from GHGRP with the corresponding ORIS code from the crosswalk. I chose an inner join because I only wanted facilities that appear in both datasets. The resulting dataframe, ghg_with_oris, now has the emissions data and ORIS codes side by side, which will let me link these facilities to electricity generation data from the EIA.

Finally, I inspected the merged dataset:

Checking .shape gave me the number of rows and columns, so I could see how many facilities were successfully matched.

Using .head() allowed me to preview the first few rows and make sure the merge worked as expected. Now I have a combined dataset where each GHGRP facility’s emissions are linked to its ORIS code, ready for merging with the EIA power generation data.

ghg_with_oris = pd.merge(
    ghg_small, 
    cross_small, 
    on="Facility_Id", 
    how="inner"
)

print("Merged GHGRP:", ghg_with_oris.shape)
ghg_with_oris.head()

Merged GHGRP: (1643, 5)

	Facility_Id	Facility Name	State	Total reported direct emissions	ORIS CODE
0	1000112.0	23rd and 3rd	NY	31916.132	7910
1	1001106.0	48th Street Peaking Station	MI	39263.556	7258
2	1000766.0	59th Street	NY	174671.338	2503
3	1000767.0	74th Street	NY	251175.596	2504
4	1001033.0	A B Brown Generating Station	IN	1982200.776	6137

Next, I loaded the EIA Form 923 dataset, which contains monthly electricity generation data for U.S. power plants. I used pd.read_csv() to read the CSV file and this time, I skipped the first five rows because they contained extra header information rather than actual column names.

I then inspected the columns with eia.columns to understand what variables were available.

eia = pd.read_csv(
    "/Users/ariakovalovich/Documents/SIPA MIA Coursework/aak2257.github.io/"
    "Project 2 Datasets/EIA923_Schedules_2_3_4_5_M_12_2023_Final_Revision.csv",
    skiprows=5)
eia.columns

/var/folders/pl/4_5qmb416_v7ywlr51mndq9h0000gn/T/ipykernel_41382/923542817.py:1: DtypeWarning:

Columns (5) have mixed types. Specify dtype option on import or set low_memory=False.

Index(['Plant Id', 'Combined Heat And\nPower Plant', 'Nuclear Unit Id',
       'Plant Name', 'Operator Name', 'Operator Id', 'Plant State',
       'Census Region', 'NERC Region', 'Reserved', 'NAICS Code',
       'EIA Sector Number', 'Sector Name', 'Reported\nPrime Mover',
       'Reported\nFuel Type Code', 'MER\nFuel Type Code',
       'Balancing\nAuthority Code', 'Respondent\nFrequency',
       'Physical\nUnit Label', 'Quantity\nJanuary', 'Quantity\nFebruary',
       'Quantity\nMarch', 'Quantity\nApril', 'Quantity\nMay', 'Quantity\nJune',
       'Quantity\nJuly', 'Quantity\nAugust', 'Quantity\nSeptember',
       'Quantity\nOctober', 'Quantity\nNovember', 'Quantity\nDecember',
       'Elec_Quantity\nJanuary', 'Elec_Quantity\nFebruary',
       'Elec_Quantity\nMarch', 'Elec_Quantity\nApril', 'Elec_Quantity\nMay',
       'Elec_Quantity\nJune', 'Elec_Quantity\nJuly', 'Elec_Quantity\nAugust',
       'Elec_Quantity\nSeptember', 'Elec_Quantity\nOctober',
       'Elec_Quantity\nNovember', 'Elec_Quantity\nDecember',
       'MMBtuPer_Unit\nJanuary', 'MMBtuPer_Unit\nFebruary',
       'MMBtuPer_Unit\nMarch', 'MMBtuPer_Unit\nApril', 'MMBtuPer_Unit\nMay',
       'MMBtuPer_Unit\nJune', 'MMBtuPer_Unit\nJuly', 'MMBtuPer_Unit\nAugust',
       'MMBtuPer_Unit\nSeptember', 'MMBtuPer_Unit\nOctober',
       'MMBtuPer_Unit\nNovember', 'MMBtuPer_Unit\nDecember',
       'Tot_MMBtu\nJanuary', 'Tot_MMBtu\nFebruary', 'Tot_MMBtu\nMarch',
       'Tot_MMBtu\nApril', 'Tot_MMBtu\nMay', 'Tot_MMBtu\nJune',
       'Tot_MMBtu\nJuly', 'Tot_MMBtu\nAugust', 'Tot_MMBtu\nSeptember',
       'Tot_MMBtu\nOctober', 'Tot_MMBtu\nNovember', 'Tot_MMBtu\nDecember',
       'Elec_MMBtu\nJanuary', 'Elec_MMBtu\nFebruary', 'Elec_MMBtu\nMarch',
       'Elec_MMBtu\nApril', 'Elec_MMBtu\nMay', 'Elec_MMBtu\nJune',
       'Elec_MMBtu\nJuly', 'Elec_MMBtu\nAugust', 'Elec_MMBtu\nSeptember',
       'Elec_MMBtu\nOctober', 'Elec_MMBtu\nNovember', 'Elec_MMBtu\nDecember',
       'Netgen\nJanuary', 'Netgen\nFebruary', 'Netgen\nMarch', 'Netgen\nApril',
       'Netgen\nMay', 'Netgen\nJune', 'Netgen\nJuly', 'Netgen\nAugust',
       'Netgen\nSeptember', 'Netgen\nOctober', 'Netgen\nNovember',
       'Netgen\nDecember', 'Total Fuel Consumption\nQuantity',
       'Electric Fuel Consumption\nQuantity', 'Total Fuel Consumption\nMMBtu',
       'Elec Fuel Consumption\nMMBtu', 'Net Generation\n(Megawatthours)',
       'YEAR'],
      dtype='object')

From the full dataset, I selected only the columns relevant for my analysis:

Plant Id → the unique identifier for each power plant (ORIS code)
Plant Name → the name of the power plant
Plant State → the state where the plant is located
Net Generation (MWh) → the total electricity generated by the plant
Reported Fuel Type Code → the primary fuel used for generation

I stored these columns in a smaller DataFrame called eia_small to make subsequent analysis and merging easier:

To make the column names consistent and easier to work with, I renamed “Plant Id” to “ORIS CODE.”

Because some plants reported multiple rows of generation data (because some power plants contain multiple generating units), I aggregated the data by ORIS Code to get a single row per plant.

I summed the total generation and kept the first value for plant name, state, and fuel type.

Finally, I inspected the aggregated dataset with .head() to make sure everything looked correct.

eia_small = eia[
    [
        "Plant Id",
        "Plant Name",
        "Plant State",
        "Net Generation\n(Megawatthours)",
        "Reported\nFuel Type Code",
    ]
]
eia_small = eia_small.rename(columns={"Plant Id": "ORIS CODE"})
eia_small.columns
eia_small.head()    

eia_agg = eia_small.groupby("ORIS CODE", as_index=False).agg(
    {
        "Net Generation\n(Megawatthours)": "sum",
        "Plant Name": "first",
        "Plant State": "first",
        "Reported\nFuel Type Code": "first",
    }
)

eia_agg.head()

	ORIS CODE	Net Generation\n(Megawatthours)	Plant Name	Plant State	Reported\nFuel Type Code
0	1	3192,74851183	Sand Point	AK	DFO
1	2	138,682	Bankhead Dam	AL	WAT
2	3	2,224,7844,324,0561,529,006597,9500	Barry	AL	NG
3	4	514,504	Walter Bouldin Dam	AL	WAT
4	7	0	Gadsden	AL	NG

I tried merging the GHGRP emissions data with the EIA generation data, but was getting an error message that the ORIS Code column in ghg_with_oris was a different type (dtype(‘O’) than the EIA ORIS column, which was type int. I checked the dataset and saw that some rows of the ORIS column contained “No Match”.

So I converted the Oris columns in both datasets to numbers using pd.to_numeric(), specifying errors=”coerce” to turn any invalid or missing entries into NaN.

ghg_with_oris.columns
ghg_with_oris 
ghg_with_oris["ORIS CODE"].dtype
ghg_with_oris["ORIS CODE"] = pd.to_numeric(ghg_with_oris["ORIS CODE"], 
                                           errors="coerce")
eia_agg["ORIS CODE"] = pd.to_numeric(eia_agg["ORIS CODE"], errors="coerce")

ghg_with_oris.head()
eia_agg.head()

	ORIS CODE	Net Generation\n(Megawatthours)	Plant Name	Plant State	Reported\nFuel Type Code
0	1	3192,74851183	Sand Point	AK	DFO
1	2	138,682	Bankhead Dam	AL	WAT
2	3	2,224,7844,324,0561,529,006597,9500	Barry	AL	NG
3	4	514,504	Walter Bouldin Dam	AL	WAT
4	7	0	Gadsden	AL	NG

I checked types again and found that the GHG dataset was float and EIA dataset was int.

print("GHG ORIS dtype:", ghg_with_oris["ORIS CODE"].dtype)          
print("EIA ORIS dtype:", eia_agg["ORIS CODE"].dtype)

GHG ORIS dtype: float64
EIA ORIS dtype: int64

Next, I addressed the issue of missing or invalid ORIS codes. In the EIA dataset, some rows converted to NaN when I used pd.to_numeric() because they contained non-numeric values (like “No Match”). To clean this up, I dropped any rows where “ORIS CODE” was missing.

Finally, I converted the remaining ORIS codes to integers to ensure a consistent format for merging.

eia_agg = eia_agg.dropna(subset=["ORIS CODE"])
eia_agg["ORIS CODE"] = eia_agg["ORIS CODE"].astype(int)
ghg_with_oris = ghg_with_oris.dropna(subset=["ORIS CODE"])
ghg_with_oris["ORIS CODE"] = ghg_with_oris["ORIS CODE"].astype(int)


print("GHG ORIS dtype:", ghg_with_oris["ORIS CODE"].dtype)
print("EIA ORIS dtype:", eia_agg["ORIS CODE"].dtype)

GHG ORIS dtype: int64
EIA ORIS dtype: int64

With the ORIS codes cleaned and formatted consistently in both datasets, I was ready to merge the GHGRP emissions data with the EIA generation data. I performed an inner merge on the “ORIS CODE” column. By merging on “ORIS CODE”, I was able to link each facility’s emissions to its electricity generation and fuel type.

I chose an inner join because I only wanted to include power plants that appear in both datasets, ensuring that every row in the resulting DataFrame has both emissions and generation data.

After merging, I checked the shape of the resulting DataFrame. This told me how many power plants were successfully matched. Finally, I used .head() to inspect the first few rows of final_merge and confirm that the merge worked as expected.

final_merge = pd.merge(ghg_with_oris, eia_agg, on="ORIS CODE", how="inner")

print("Final merged shape:", final_merge.shape)
final_merge.head()

Final merged shape: (1493, 9)

	Facility_Id	Facility Name	State	Total reported direct emissions	ORIS CODE	Net Generation\n(Megawatthours)	Plant Name	Plant State	Reported\nFuel Type Code
0	1000112.0	23rd and 3rd	NY	31916.132	7910	52,540	Joseph J Seymour Power Project	NY	NG
1	1000766.0	59th Street	NY	174671.338	2503	0612	59th Street	NY	KER
2	1000767.0	74th Street	NY	251175.596	2504	176	74th Street	NY	KER
3	1001033.0	A B Brown Generating Station	IN	1982200.776	6137	41231,3631,719,05932,335	A B Brown	IN	DFO
4	1013527.0	A.J. Mihm Generating Station	MI	106545.348	61391	207,891	A.J. Mihm Generating Station	MI	NG

I’ve worked in energy and climate policy for almost ten years, and even I was confused by the EIA’s fuel type classification system. I looked up a key online, and made it into a Python dictionary, grouping the fuel types into broader categories that a layperson should be able to understand, to make the eventual visualization by fuel type easier to interpret.

Coal: BIT, LIG, SUB, WC, RC

Oil: DFO, RFO, KER, JF, PC, WO

Gas: NG, BFG, COG, LFG, OG, PRG, SGC

Biomass: AB, BLQ, LFG, MSW, WDL, WDS

Renewables: GEO, H2, NUC, SUN, WND

Other: OTH

I used the .map() function to create a new column “Fuel_Category” in my merged dataset, final_merge, where each plant’s fuel code is replaced with its broader category, which I planned to use for coloring each point on the scatter plot by broad fuel type to visualize emissions across different types of power generation.

fuel_map = {
    # Coal
    "BIT": "Coal",
    "LIG": "Coal",
    "SUB": "Coal",
    "WC": "Coal",
    "RC": "Coal",
    # Oil
    "DFO": "Oil",
    "RFO": "Oil",
    "KER": "Oil",
    "JF": "Oil",
    "PC": "Oil",
    "WO": "Oil",
    # Gas
    "NG": "Gas",
    "BFG": "Gas",
    "COG": "Gas",
    "LFG": "Gas",
    "OG": "Gas",
    "PRG": "Gas",
    "SGC": "Gas",
    # Biomass
    "AB": "Biomass",
    "BLQ": "Biomass",
    "LFG": "Biomass",
    "MSW": "Biomass",
    "WDL": "Biomass",
    "WDS": "Biomass",
    # Renewables
    "GEO": "Clean",
    "H2": "Clean",
    "NUC": "Clean",
    "SUN": "Clean",
    "WND": "Clean",
    # Other
    "OTH": "Other",
}
final_merge["Fuel_Category"] = final_merge["Reported\nFuel Type Code"].map(fuel_map)

Now that I’ve merged the dataset and categorized the fuel types, I created a scatter plot using Plotly Express to explore the relationship between electricity generation and greenhouse gas emissions across power plants.

I plotted net generation (MWh) on the x-axis and total reported direct emissions (metric tons) on the y-axis. I used the “Fuel_Category” column to color-code the points, so it’s easy to see how emissions differ by fuel type.

The hover_data feature shows additional details for each plant when I hover over a point, including plant names, state, and the original fuel code. The labels dictionary makes the axis titles and color legend more readable.

And, to better visualize plants across the wide range of generation and emissions values, I set both axes to a logarithmic scale. Finally, I displayed the plot with fig.show().

fig = px.scatter(
    final_merge,
    x="Net Generation\n(Megawatthours)",
    y="Total reported direct emissions",
    color="Fuel_Category",  # updated categories
    hover_data=[
        "Plant Name",
        "Facility Name",
        "Plant State",
        "State",
        "Reported\nFuel Type Code",
    ],
    labels={
        "Net Generation\n(Megawatthours)": "Net Generation (MWh)",
        "Total reported direct emissions": "Total Direct Emissions (Metric Tons)",
        "Fuel_Category": "Fuel Type Category",
    },
    title="Power Plant Emissions vs. Power Generation by Fuel Type Category",
)

fig.update_layout(xaxis_type="log", yaxis_type="log")
fig.show()

Looking at the scatter plot, several patterns immediately stand out. First, as expected, there is a positive correlation between net generation and total emissions—plants that generate more electricity generally emit more greenhouse gases. However, the relationship is not perfectly linear, which reflects differences in plant efficiency and the emissions intensities of different fuels.

It is easy to see that coal, gas, and oil plants have the highest emissions, and this seems to hold through for a wide range of output generation, which I didn’t expect.

I would have expected the relactionships to be more linear within a fuel type, wtih different clusters according to the relative emissions intensities of each fuel. For example, I expected coal and oil plants to dominate a high-emission, high-generation quadrant, and gas-fired plants to cluster lower in emissions for comparable generation, reflecting their lower carbon intensity relative to coal.

Instead, coal, oil, and gas plants are more scattered, with plants with the same generation exhibiting widely different emissions profiles. This could be due to the way I lumped multiple subtypes of fuel by overall type, and some of the broad categories reflect fuels of varying qualities and carbon intensities.

It’s easy to see from the graph that biomass plants tend to generate less electricity overall than the other types, and that clean plants like renewables and nuclear barely register on the graph because they produce minimal direct emissions regardless of generation.

Because I used a logarithmic scale on both axes, the plot highlights differences across plants of vastly different sizes, allowing smaller plants to be visible rather than being squashed near the origin. This makes it easier to compare emissions efficiency across fuel types, rather than just looking at the largest emitters.

In sum, while total emissions are largely driven by plant size and fuel type, there is significant variability even within categories, underscoring the importance of examining emissions intensity when evaluating the carbon efficiency of electricity generation.