Twin y axes plot with monthly average as x-axis over multiple years with xarray

Question

Thank you for taking interest in this post.

I am hoping to create a twin y axes plot illustrating 2 variables (eg. y1=snowmelt and y2=discharge run-off). Due to my very limited knowledge in python as I have just started to learn coding, I am very confused as to how I can go about it.

I am able to construct a twin y axes plot and a plot representing monthly average of multiple years. However, I am not sure how to combine these matplotlib codes together and create a twin y axes plot with monthly average as x-axis over multiple years.

Code for monthly average

disda = xr.open_mfdataset(sorted(['1980.nc','1981.nc', '1982.nc', '1983.nc','1984.nc','1985.nc','1986.nc','1987.nc', '1988.nc','1989.nc', '1990.nc', '1991.nc', '1992.nc', '1993.nc', '1994.nc', '1996.nc', '1997.nc', '1998.nc','1999.nc']))
snowdepth_we = xr.open_dataarray('Snow depth water equivalent.nc')


disyrs = disda.sel(time=slice('1981','1999'))
snyrs = snowdepth_we.sel(time=slice('1981','1999'))


dismonthlymean = disyrs.dis24

dislatlonmean = dismonthlymean.mean(dim=['latitude','longitude']).resample(time='M').sum()
snowlatlonmean = snmonthlymean.mean(dim=['latitude','longitude'])


disgroupby = dislatlonmean.groupby("time.month").mean("time")
sngroupby = snowlatlonmean.groupby("time.month").mean("time")

#graph commands
myfig, myax = plt.subplots(figsize=(12,6))

    
TYs = np.unique(dislatlonmean["time.year"])
disgroupby.plot.line('b-', color='blue', linestyle='-', linewidth=4, label='Discharge Mean')
for YY in TYs: 
    plt.plot(np.arange(1,13), dislatlonmean.sel(time=YY.astype("str")), 'b-', color='blue', alpha=0.2)        

TYs = np.unique(dislatlonmean["time.year"])
sngroupby.plot.line('b-', color='red', linestyle='-', linewidth=4, label='Snowdepth Mean')
for YY in TYs: 
    plt.plot(np.arange(1,13), dislatlonmean.sel(time=YY.astype("str")), 'b-', color='blue', alpha=0.2)        
        

    
myax.set_title('Western Himalayas Snow Depth and River Indus Discharge Run-off 1981-1999')
myax.set_ylabel('m of water equivalent')
myax.set_xlabel('Month')
myax.set_xticks(range(1, 13))
myax.set_xticklabels(['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sept','Oct','Nov','Dec'])
myax.grid()
myax.legend(['Discharge Mean','Snowdepth'])

since this is not a twin y-axes plot, i was unable to plot 2 variables with different measurements together where one variable became a flat line in the graph, underneath the other graph

on the other hand, i was able to create a twin y-axes plot using these commands

disda_dis248199 = disda_dis24.loc['1981':'1999'].mean(dim=['latitude','longitude']).resample(time='M').mean()
snow8199 = snowdepth_we.loc['1981':'1999'].mean(dim=['latitude','longitude']).resample(time='M').mean()

x = disda_dis248199.time
y1 = snow8199
y2 = disda_dis248199

#plotting commands
fig, ax1 = plt.subplots(figsize=(14,8))

ax1.set_xlabel('Year')
ax1.set_ylabel('Snow depth water equivalent (m of w.e)')
ax1.plot(x, y1, color='red', label='River Discharge')
ax1.tick_params(axis='y', labelcolor='red')

ax2 = ax1.twinx()  # instantiate a second axes that shares the same x-axis

ax2.set_ylabel('River Discharge (m3s−1)')  # we already handled the x-label with ax1
ax2.plot(x, y2, color='blue',alpha=0.4, label='Snow depth water equivalent')
ax2.tick_params(axis='y', labelcolor='blue')

fig.legend(loc="upper right")
ax1.set_title('Western Himalayas Snow Depth and River Indus Tributaries Comparison 1981 - 2016')
ax2.grid()

I have tried to merge these 2 commands together, but to no avail.

As a result, I am wondering if there's any way I can plot a twin y-axes plot like the 2nd graph but with months as the x-axis and both lines on the graph will be representing 2 different variables.

In addition, as you can see in graph 1, somehow the legend only shows the discharge mean (bold lines) and actual yearly averages from 1981-1999 (thin lines), therefore, I am wondering how would I be able to fix this?

I greatly appreciate your help!! Please let me know if you require any more information, I will respond as soon as possible.

Answer 1

Since you only gave part of your code, I cannot run it and thus also cannot change it for you. However, I will try to give an explanation on what should be changed in the first piece of code.

As a first step, you want to define a second Axes object similarly to how you did it in the second piece of code:

myax2 = myax.twinx()

Since you now have two Axes instances, you'll need to specify on which Axes instance (either myax or myax2) you want to plot the data. This means that instead of

disgroupby.plot.line('b-', color='blue', linestyle='-', linewidth=4, label='Discharge Mean')
sngroupby.plot.line('b-', color='red', linestyle='-', linewidth=4, label='Snowdepth Mean')

you should do something such as

disgroupby.plot.line('b-', ax=myax, color='blue', linestyle='-', linewidth=4, label='Discharge Mean')
sngroupby.plot.line('b-', ax=myax2, color='red', linestyle='-', linewidth=4, label='Snowdepth Mean')

Note the addition of ax=... in the list of (keyword) arguments. See this part of xarray's documentation for an example. In addition, you no longer want to use plt.plot, as this does not specify on which Axes instance to plot the data. Instead, you want to call the plot method of the Axes instance on which you want to plot the data. As an example, the first occurence of

plt.plot(np.arange(1,13), dislatlonmean.sel(time=YY.astype("str")), 'b-', color='blue', alpha=0.2)  

should be changed to

myax.plot(np.arange(1,13), dislatlonmean.sel(time=YY.astype("str")), 'b-', color='blue', alpha=0.2)

and the second occurence should be changed to

myax2.plot(np.arange(1,13), dislatlonmean.sel(time=YY.astype("str")), 'b-', color='blue', alpha=0.2)

Setting labels, titles, etc. is similar to how you did it in the second piece of code.

Small note about your plt.plot statements, you have ..., 'b-', color='blue', ... in the list of arguments. While this works, you specify the color twice. You could either remove the color='blue' part or change 'b-' to '-'.