Historical time-series analysis of Many Glacier, Montana¶

This analysis looks at climate data from Many Glacier, MT using the NOAA/NCEI API.

In [1]:
# Import python libraries

import pandas as pd  # working with datafames 
import geopandas as gpd  # working with vector data
import osmnx as ox  # creating maps with open street maps

import holoviews as hv  # creating interactive plots
import hvplot.pandas

import numpy as np  # for adding trendline to plot
import matplotlib.pyplot as plt  # for plotting
In [2]:
# Build endpoint request from NOAA/NCEI for Many Glacier, MT

many_glacier_url = (
    "https://www.ncei.noaa.gov/access/services/data/v1"
    "?dataset=daily-summaries"
    "&stations=USS0013A27S"
    "&startDate=1976-09-30"
    "&endDate=2024-04-25"
    "&dataTypes=TOBS,PRCP,WESD"
    "&includeStationName=true"
    "&includeStationLocation=1"
    "&units=standard")

many_glacier_url
Out[2]:
'https://www.ncei.noaa.gov/access/services/data/v1?dataset=daily-summaries&stations=USS0013A27S&startDate=1976-09-30&endDate=2024-04-25&dataTypes=TOBS,PRCP,WESD&includeStationName=true&includeStationLocation=1&units=standard'
In [3]:
# Open up data using pandas
many_glacier_df = pd.read_csv(many_glacier_url,
                              parse_dates=True,
                              index_col='DATE',
                              na_values=['NaN'])

# Clean up data (drop columns that aren't varying with time)
many_glacier_df = many_glacier_df[['PRCP','TOBS','WESD']]
many_glacier_df
Out[3]:
PRCP TOBS WESD
DATE
1976-09-30 NaN NaN 0.0
1976-10-01 NaN NaN 0.0
1976-10-02 NaN NaN 0.0
1976-10-03 NaN NaN 0.0
1976-10-04 NaN NaN 0.0
... ... ... ...
2024-04-21 0.0 37.0 4.2
2024-04-22 0.0 30.0 3.5
2024-04-23 0.0 43.0 2.5
2024-04-24 0.0 45.0 1.4
2024-04-25 0.0 36.0 1.0

17375 rows × 3 columns

In [4]:
# Plot WESD data using .plot()
many_glacier_df.plot(y='WESD')
Out[4]:
<Axes: xlabel='DATE'>
No description has been provided for this image
In [5]:
# Plot temperature data using .hvplot()
many_glacier_df.hvplot(y='TOBS')
Out[5]:
In [6]:
# Resample to look annual mean 
many_glacier_annual_avg = many_glacier_df.resample('YE').mean()
many_glacier_annual_avg
Out[6]:
PRCP TOBS WESD
DATE
1976-12-31 NaN NaN 1.236559
1977-12-31 NaN NaN 3.251233
1978-12-31 0.073913 16.958333 6.226301
1979-12-31 0.100274 30.415473 5.050137
1980-12-31 0.157377 30.069364 5.042077
1981-12-31 0.110137 29.729008 3.637260
1982-12-31 0.150137 26.740741 6.632603
1983-12-31 0.108767 28.765258 4.653151
1984-12-31 0.138525 30.477912 3.904098
1985-12-31 0.132329 28.224658 5.284110
1986-12-31 0.150411 33.024658 2.811781
1987-12-31 0.106027 34.741047 3.114795
1988-12-31 0.125137 33.284153 2.660109
1989-12-31 0.171507 31.958217 4.486301
1990-12-31 0.186849 31.986301 4.917260
1991-12-31 0.140548 33.183562 7.052329
1992-12-31 0.117760 33.980822 2.566120
1993-12-31 0.106027 31.541547 3.400274
1994-12-31 0.127123 34.595989 4.601096
1995-12-31 0.184932 33.128767 2.691233
1996-12-31 0.156557 30.441096 4.346175
1997-12-31 0.133425 34.379121 7.852877
1998-12-31 0.127945 36.800000 3.651233
1999-12-31 0.138356 37.393939 5.363562
2000-12-31 0.100000 35.227397 3.880874
2001-12-31 0.098082 36.912329 3.098082
2002-12-31 0.152603 34.783562 5.455342
2003-12-31 0.110685 37.751381 2.910959
2004-12-31 0.129781 37.265027 2.619399
2005-12-31 0.122740 37.824658 1.066575
2006-12-31 0.142192 38.323288 2.756164
2007-12-31 0.113151 38.583562 2.555616
2008-12-31 0.118579 36.633880 4.435246
2009-12-31 0.116712 35.391781 3.484110
2010-12-31 0.120548 37.452055 2.923014
2011-12-31 0.149041 36.564384 6.660548
2012-12-31 0.163661 37.961749 4.326503
2013-12-31 0.119178 38.057534 3.069589
2014-12-31 0.144384 37.046575 4.480000
2015-12-31 0.118630 39.917808 1.161370
2016-12-31 0.134153 39.412568 1.610109
2017-12-31 0.129041 37.030137 3.955616
2018-12-31 0.100274 37.098630 5.406849
2019-12-31 0.112329 34.886427 2.983288
2020-12-31 0.125956 38.041436 4.150273
2021-12-31 0.116986 37.831025 3.002466
2022-12-31 0.120548 36.825485 4.346301
2023-12-31 0.103014 38.564384 3.028219
2024-12-31 0.132759 26.370690 4.487931
In [7]:
# Subset the data to look at 1980-2020
many_glacier_1980_2020 = many_glacier_annual_avg["1980":"2020"]
many_glacier_1980_2020
Out[7]:
PRCP TOBS WESD
DATE
1980-12-31 0.157377 30.069364 5.042077
1981-12-31 0.110137 29.729008 3.637260
1982-12-31 0.150137 26.740741 6.632603
1983-12-31 0.108767 28.765258 4.653151
1984-12-31 0.138525 30.477912 3.904098
1985-12-31 0.132329 28.224658 5.284110
1986-12-31 0.150411 33.024658 2.811781
1987-12-31 0.106027 34.741047 3.114795
1988-12-31 0.125137 33.284153 2.660109
1989-12-31 0.171507 31.958217 4.486301
1990-12-31 0.186849 31.986301 4.917260
1991-12-31 0.140548 33.183562 7.052329
1992-12-31 0.117760 33.980822 2.566120
1993-12-31 0.106027 31.541547 3.400274
1994-12-31 0.127123 34.595989 4.601096
1995-12-31 0.184932 33.128767 2.691233
1996-12-31 0.156557 30.441096 4.346175
1997-12-31 0.133425 34.379121 7.852877
1998-12-31 0.127945 36.800000 3.651233
1999-12-31 0.138356 37.393939 5.363562
2000-12-31 0.100000 35.227397 3.880874
2001-12-31 0.098082 36.912329 3.098082
2002-12-31 0.152603 34.783562 5.455342
2003-12-31 0.110685 37.751381 2.910959
2004-12-31 0.129781 37.265027 2.619399
2005-12-31 0.122740 37.824658 1.066575
2006-12-31 0.142192 38.323288 2.756164
2007-12-31 0.113151 38.583562 2.555616
2008-12-31 0.118579 36.633880 4.435246
2009-12-31 0.116712 35.391781 3.484110
2010-12-31 0.120548 37.452055 2.923014
2011-12-31 0.149041 36.564384 6.660548
2012-12-31 0.163661 37.961749 4.326503
2013-12-31 0.119178 38.057534 3.069589
2014-12-31 0.144384 37.046575 4.480000
2015-12-31 0.118630 39.917808 1.161370
2016-12-31 0.134153 39.412568 1.610109
2017-12-31 0.129041 37.030137 3.955616
2018-12-31 0.100274 37.098630 5.406849
2019-12-31 0.112329 34.886427 2.983288
2020-12-31 0.125956 38.041436 4.150273
In [8]:
# Resetting the index
many_glacier_1980_2020 = many_glacier_1980_2020.reset_index()
many_glacier_1980_2020
Out[8]:
DATE PRCP TOBS WESD
0 1980-12-31 0.157377 30.069364 5.042077
1 1981-12-31 0.110137 29.729008 3.637260
2 1982-12-31 0.150137 26.740741 6.632603
3 1983-12-31 0.108767 28.765258 4.653151
4 1984-12-31 0.138525 30.477912 3.904098
5 1985-12-31 0.132329 28.224658 5.284110
6 1986-12-31 0.150411 33.024658 2.811781
7 1987-12-31 0.106027 34.741047 3.114795
8 1988-12-31 0.125137 33.284153 2.660109
9 1989-12-31 0.171507 31.958217 4.486301
10 1990-12-31 0.186849 31.986301 4.917260
11 1991-12-31 0.140548 33.183562 7.052329
12 1992-12-31 0.117760 33.980822 2.566120
13 1993-12-31 0.106027 31.541547 3.400274
14 1994-12-31 0.127123 34.595989 4.601096
15 1995-12-31 0.184932 33.128767 2.691233
16 1996-12-31 0.156557 30.441096 4.346175
17 1997-12-31 0.133425 34.379121 7.852877
18 1998-12-31 0.127945 36.800000 3.651233
19 1999-12-31 0.138356 37.393939 5.363562
20 2000-12-31 0.100000 35.227397 3.880874
21 2001-12-31 0.098082 36.912329 3.098082
22 2002-12-31 0.152603 34.783562 5.455342
23 2003-12-31 0.110685 37.751381 2.910959
24 2004-12-31 0.129781 37.265027 2.619399
25 2005-12-31 0.122740 37.824658 1.066575
26 2006-12-31 0.142192 38.323288 2.756164
27 2007-12-31 0.113151 38.583562 2.555616
28 2008-12-31 0.118579 36.633880 4.435246
29 2009-12-31 0.116712 35.391781 3.484110
30 2010-12-31 0.120548 37.452055 2.923014
31 2011-12-31 0.149041 36.564384 6.660548
32 2012-12-31 0.163661 37.961749 4.326503
33 2013-12-31 0.119178 38.057534 3.069589
34 2014-12-31 0.144384 37.046575 4.480000
35 2015-12-31 0.118630 39.917808 1.161370
36 2016-12-31 0.134153 39.412568 1.610109
37 2017-12-31 0.129041 37.030137 3.955616
38 2018-12-31 0.100274 37.098630 5.406849
39 2019-12-31 0.112329 34.886427 2.983288
40 2020-12-31 0.125956 38.041436 4.150273
In [9]:
# Remove year from DATE column and add as new variable
many_glacier_1980_2020['YEAR'] = many_glacier_1980_2020['DATE'].dt.year
many_glacier_1980_2020
Out[9]:
DATE PRCP TOBS WESD YEAR
0 1980-12-31 0.157377 30.069364 5.042077 1980
1 1981-12-31 0.110137 29.729008 3.637260 1981
2 1982-12-31 0.150137 26.740741 6.632603 1982
3 1983-12-31 0.108767 28.765258 4.653151 1983
4 1984-12-31 0.138525 30.477912 3.904098 1984
5 1985-12-31 0.132329 28.224658 5.284110 1985
6 1986-12-31 0.150411 33.024658 2.811781 1986
7 1987-12-31 0.106027 34.741047 3.114795 1987
8 1988-12-31 0.125137 33.284153 2.660109 1988
9 1989-12-31 0.171507 31.958217 4.486301 1989
10 1990-12-31 0.186849 31.986301 4.917260 1990
11 1991-12-31 0.140548 33.183562 7.052329 1991
12 1992-12-31 0.117760 33.980822 2.566120 1992
13 1993-12-31 0.106027 31.541547 3.400274 1993
14 1994-12-31 0.127123 34.595989 4.601096 1994
15 1995-12-31 0.184932 33.128767 2.691233 1995
16 1996-12-31 0.156557 30.441096 4.346175 1996
17 1997-12-31 0.133425 34.379121 7.852877 1997
18 1998-12-31 0.127945 36.800000 3.651233 1998
19 1999-12-31 0.138356 37.393939 5.363562 1999
20 2000-12-31 0.100000 35.227397 3.880874 2000
21 2001-12-31 0.098082 36.912329 3.098082 2001
22 2002-12-31 0.152603 34.783562 5.455342 2002
23 2003-12-31 0.110685 37.751381 2.910959 2003
24 2004-12-31 0.129781 37.265027 2.619399 2004
25 2005-12-31 0.122740 37.824658 1.066575 2005
26 2006-12-31 0.142192 38.323288 2.756164 2006
27 2007-12-31 0.113151 38.583562 2.555616 2007
28 2008-12-31 0.118579 36.633880 4.435246 2008
29 2009-12-31 0.116712 35.391781 3.484110 2009
30 2010-12-31 0.120548 37.452055 2.923014 2010
31 2011-12-31 0.149041 36.564384 6.660548 2011
32 2012-12-31 0.163661 37.961749 4.326503 2012
33 2013-12-31 0.119178 38.057534 3.069589 2013
34 2014-12-31 0.144384 37.046575 4.480000 2014
35 2015-12-31 0.118630 39.917808 1.161370 2015
36 2016-12-31 0.134153 39.412568 1.610109 2016
37 2017-12-31 0.129041 37.030137 3.955616 2017
38 2018-12-31 0.100274 37.098630 5.406849 2018
39 2019-12-31 0.112329 34.886427 2.983288 2019
40 2020-12-31 0.125956 38.041436 4.150273 2020
In [10]:
# Plot TOBS using .plot()
many_glacier_1980_2020.plot(y='TOBS',
                            x='YEAR')
Out[10]:
<Axes: xlabel='YEAR'>
No description has been provided for this image
In [38]:
# Plot PRCP using .hvplot()
many_glacier_1980_2020.hvplot(y='PRCP',
                              x='YEAR',
                              title='Change in Precipitation at Many Glacier, MT - 1980-2020')
Out[38]:

Plotting with Matplotlib and ChatGPT¶

Read more about a matplotlib at our open textbook.

In [39]:
# Define our figure and axis objects 
fig, ax = plt.subplots(figsize=(6,4))

# Plot TOBS vs. YEAR as scatter plot
# Plot scatter set
ax.scatter(many_glacier_1980_2020['YEAR'],
           many_glacier_1980_2020['TOBS'],
           color='white',
           edgecolor='black')

# Plot line plot
ax.plot(many_glacier_1980_2020['YEAR'],
        many_glacier_1980_2020['TOBS'],
        color='grey')

# Add title and axis label
ax.set(title="Mean Annual Temperature\nMany Glacier, MT (1980-2020)",
       ylabel="Temperature (F)")
Out[39]:
[Text(0.5, 1.0, 'Mean Annual Temperature\nMany Glacier, MT (1980-2020)'),
 Text(0, 0.5, 'Temperature (F)')]
No description has been provided for this image
In [40]:
# From ChatGPT

# Define our figure and axis objects 
fig, ax = plt.subplots(figsize=(6,4))

# Compute linear regression
x = many_glacier_1980_2020['YEAR']
y = many_glacier_1980_2020['TOBS']

# Compute the slope (m) and intercept (b) of the line y = mx + b
m, b = np.polyfit(x, y, 1)

# Plot TOBS vs. YEAR as scatter plot
ax.scatter(x, y, color='white', edgecolor='black')

# Plot trend line
ax.plot(x, m*x + b, color='blue', label=f'Trend Line (R-squared = {np.corrcoef(x, y)[0,1]**2:.2f})')

# Add legend
ax.legend()

# Add title and axis label
ax.set(title="Mean Annual Temperature\nMany Glacier, MT (1980-2020)",
       ylabel="Temperature (F)")
Out[40]:
[Text(0.5, 1.0, 'Mean Annual Temperature\nMany Glacier, MT (1980-2020)'),
 Text(0, 0.5, 'Temperature (F)')]
No description has been provided for this image
In [41]:
# Define our figure and axis objects 
fig, ax = plt.subplots(figsize=(6,4))

# Compute linear regression
x = many_glacier_1980_2020['YEAR']
y = many_glacier_1980_2020['TOBS']

# Compute the slope (m) and intercept (b) of the line y = mx + b
m, b = np.polyfit(x, y, 1)

# Plot TOBS vs. YEAR as scatter plot
ax.scatter(x, y, color='white', edgecolor='black')

# Plot trend line
ax.plot(x, m*x + b, color='blue', label=f'Trend Line (Slope = {m:.2f})')

# Add legend
ax.legend()

# Add title and axis label
ax.set(title="Mean Annual Temperature\nMany Glacier, MT (1980-2020)",
       ylabel="Temperature (F)")
Out[41]:
[Text(0.5, 1.0, 'Mean Annual Temperature\nMany Glacier, MT (1980-2020)'),
 Text(0, 0.5, 'Temperature (F)')]
No description has been provided for this image
In [42]:
# Define our figure and axis objects 
fig, ax = plt.subplots(figsize=(6,4))

# Compute linear regression
x = many_glacier_1980_2020['YEAR']
y = many_glacier_1980_2020['TOBS']

# Compute the slope (m) and intercept (b) of the line y = mx + b
m, b = np.polyfit(x, y, 1)

# Define degree symbol
degree_symbol = '\u00B0'

# Plot TOBS vs. YEAR as scatter plot
ax.scatter(x, y, color='white', edgecolor='black')

# Plot line graph
ax.plot(many_glacier_1980_2020['YEAR'],
        many_glacier_1980_2020['TOBS'],
        color='grey')

# Plot trend line
ax.plot(x, m*x + b, color='tomato', label=f'Trend Line (Slope = {m:.2f} {degree_symbol}F/yr)')

# Add legend
ax.legend()

# Add title and axis label
ax.set(title="Mean Annual Temperature\nMany Glacier, MT (1980-2020)",
       ylabel="Temperature (\u00B0F)")

plt.show()
No description has been provided for this image
In [11]:
# Define our figure and axis objects 
fig, ax = plt.subplots(figsize=(6,4))

# Compute linear regression
x = many_glacier_1980_2020['YEAR']
y = many_glacier_1980_2020['WESD']

# Compute the slope (m) and intercept (b) of the line y = mx + b
m, b = np.polyfit(x, y, 1)

# Define degree symbol
degree_symbol = '\u00B0'

# Plot TOBS vs. YEAR as scatter plot
ax.scatter(x, y, color='white', edgecolor='black')

# Plot line graph
ax.plot(many_glacier_1980_2020['YEAR'],
        many_glacier_1980_2020['WESD'],
        color='grey')

# Plot trend line
ax.plot(x, m*x + b, color='tomato', label=f'Trend Line (Slope = {m:.2f} {degree_symbol}F/yr)')

# Add legend
ax.legend()

# Add title and axis label
ax.set(title="Mean Annual Snow Water Equivalence\nMany Glacier, MT (1980-2020)",
       ylabel="Water Equivalence (in.)")

plt.show()
No description has been provided for this image

Temperatures at Many Glacier are increasing.¶

Over the last 40 years (1980-2020) the average annual temperature has increased at a rate of 0.23 deg Fahrenheit per year.

Adding a map of Glacier National Park using Open Street Maps and ChatGPT¶

In [43]:
# Define the place name (in this case, Glacier National Park)
place_name = "Glacier National Park"

# Download the boundary shapefile
gla_boundary = ox.geocode_to_gdf(place_name)
gla_boundary
Out[43]:
geometry bbox_north bbox_south bbox_east bbox_west place_id osm_type osm_id lat lon class type place_rank importance addresstype name display_name
0 POLYGON ((-114.47550 49.00109, -114.47548 49.0... 49.001109 48.233699 -113.241617 -114.475501 31377934 relation 1242641 48.617486 -113.760885 boundary national_park 25 0.514759 national_park Glacier National Park Glacier National Park, Montana, 59916, United ...
In [44]:
# Map using .explore()
gla_boundary.explore()
Out[44]:
Make this Notebook Trusted to load map: File -> Trust Notebook
In [45]:
# Plot Glacier National Park boundary
gla_map = gla_boundary.hvplot(
    # Give the map a descriptive title
    title="Glacier National Park, Montana, USA",
    # Add a basemap
    geo=True, tiles='EsriImagery',
    # Change the colors
    fill_color='white', fill_alpha=0.2,
    line_color='skyblue', line_width=5,
    # Change the image size
    frame_width=400, frame_height=400)

# Save the map as a file to put on the web
hv.save(gla_map, 'glacier-national-park.html')

# Display the map
gla_map

WARNING:bokeh.core.validation.check:W-1005 (FIXED_SIZING_MODE): 'fixed' sizing mode requires width and height to be set: figure(id='p1494', ...)
Out[45]:

Create shareable HTML of the notebook¶

In [46]:
%%capture
%%bash
jupyter nbconvert many-glacier.ipynb --to html