Geology Time Scale in Matplotlib
One of the most common things that I do as an earth scientist / geologist / stratigrapher is plotting data alongside a geologic time scale.
The idea is to show a chronostratigraphic framework of trends / data that we would want to highlight. It can be geochemical data, paleobiology data (e.g., fossil occurences, cladograms, etc), stratigraphic thickness / stacking patterns, etc.
When I first started grad school in 2017, I used an R package called deeptime. It’s a pretty nifty package, very easy to use. If you work with R, you should give this tool a try. I couldn’t recommend this enough.
However, since I converted into a devout pythonista in 2020, I had to find a way to plot geologic time scale using python. After spending days, learning about “artists” object in matplotlib, I managed to come up with a “working” (i.e., not super elegant but does its job) “go to” script that I always use whenever I need to put a geologic time scale of my data.
This is an example of the plot I made.
Here is the full code.
| """ A python Script to add geologic time scale | |
| on a matplotlib plot. | |
| """ | |
| # import libraries | |
| import numpy as np | |
| import matplotlib.pyplot as plt | |
| import pandas as pd | |
| from matplotlib.patches import Rectangle | |
| plt.rcParams['axes.linewidth'] = 1.5 | |
| # Dummy data | |
| time = np.linspace(500, 10, 200) | |
| error = np.random.normal(0, 0.75, 200) | |
| sin = np.sin(np.linspace(-3*np.pi, 3*np.pi, 200))*2 | |
| data = sin + error | |
| # Geologic Period | |
| # Boundary and colors are based on GTS 2020 and ICS | |
| gts = pd.DataFrame({ | |
| 'epoch' : ['Paleozoic', 'Mesozoic', 'Cenozoic'], | |
| 'lower_boundary' : [510, 251.9, 66], | |
| 'upper_boundary' : [251.9, 66, 0], | |
| 'color' : ['#99C08D', '#67C5CA', '#F2F91D'] | |
| }) | |
| period = pd.DataFrame({ | |
| 'epoch' : ['Cam', 'Ord', 'Sil', 'Dev', 'Car', 'Per', | |
| 'Tri', 'Jur', 'Cre', 'Pal', 'Neo', ''], | |
| 'lower_boundary' : [510, 486.9, 443.1, 419, 359.3, 298.9, 251.9, | |
| 201.4, 143.1, 66, 23.04, 2.58], | |
| 'upper_boundary' : [486.9, 443.1, 419, 359.3, 298.9, 251.9, 201.4, | |
| 143.1, 66, 23.04, 2.58, 0], | |
| 'color' : ['#7FA056', '#009270', '#B3E1B6', '#CB8C37', '#67A599', | |
| '#F04028', '#812B92', '#34B2C9', '#7FC64E', '#FD9A52', | |
| '#FFE619', '#F9F97F'] | |
| }) | |
| # Initialize figure | |
| fig = plt.figure(figsize=(8.5,5)) | |
| # Width and height ratios for gridspec | |
| w = [1] | |
| h = [1, 0.1] | |
| # Gridspec | |
| gs = fig.add_gridspec(ncols=1, nrows=2, width_ratios=w, height_ratios=h) | |
| # Initiate axes | |
| ax0 = fig.add_subplot(gs[1,0]) | |
| ax0.set_xlim(0, 510) | |
| ax = fig.add_subplot(gs[0,0], sharex=ax0) | |
| # Plot Phanerozoic Era | |
| trans = ax0.get_xaxis_transform() | |
| gts['h'] = gts['lower_boundary'] - gts['upper_boundary'] | |
| period['h'] = period['lower_boundary'] - period['upper_boundary'] | |
| rot = [90, 90, 90] | |
| for ind in gts.index: | |
| ax0.add_patch(Rectangle(xy=[gts['upper_boundary'][ind], 0], | |
| height=0.5, | |
| width=gts['h'][ind], | |
| transform=trans, ec='black', | |
| fc=gts['color'][ind])) | |
| ax0.text(gts['lower_boundary'][ind]-0.5*gts['h'][ind], 0.25, | |
| str(gts['epoch'][ind]), | |
| horizontalalignment='center', verticalalignment='center', | |
| transform = trans, fontsize=12) | |
| for ind in period.index: | |
| ax0.add_patch(Rectangle(xy=[period['upper_boundary'][ind], 0.5], | |
| height=0.5, | |
| width=period['h'][ind], | |
| transform=trans, ec='black', | |
| fc=period['color'][ind])) | |
| ax0.text(period['lower_boundary'][ind]-0.5*period['h'][ind], 0.75, | |
| str(period['epoch'][ind]), | |
| horizontalalignment='center', verticalalignment='center', | |
| transform = trans, fontsize=10) | |
| # Plot data | |
| ax.scatter(time, data, ec='k', fc='orange', s=100, | |
| label='Dummy data') | |
| # Labels and legend | |
| ax.set_ylabel('$\delta^{614}XOXO_{carb}$', fontsize=12) | |
| ax.legend() | |
| ax.tick_params(axis='x', bottom=False, labelbottom=False) | |
| ax0.tick_params(axis='y', labelleft=False, left=False) | |
| ax0.set_xlabel('Age (Ma)', fontsize=12) | |
| plt.tight_layout() | |
| plt.show() |
At some point I will post a companion post to explain each line in the code (I am not sure if it’s really necessary), but let me know if you have any questions.
Comments