Here's an explanation of what happened and some suggestions for fixing the issue.

Problem 1: Empty Plot
When you generate a panel of two subplots using ax1 and ax2, the code expects these axes objects to be created beforehand. In your case, since you didn't create any axis objects (ax1 and ax2), the plot is empty.

To fix this, simply add the following lines before creating your scatter plots:

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 5))

This creates a figure with two subplots (ax1 and ax2) arranged horizontally. The figsize parameter sets the width of each subplot.

Problem 2: Y-axis Scales
To fix this issue, you can set the same y-axis limits for both plots using the set_ylim() method. This ensures that the y-axes have the same scale and range.

ax1.set_ylim([-20, 25]) # adjust these values based on your data
ax2.set_ylim([-20, 25])

Alternatively, you can use the sharey argument when creating the subplots. This will make both y-axes have the same scale and range:

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 5), sharey=True)

Additional Suggestions

1. Consider using a single plot with two subplots (subplots(1, 2) instead of 1) to keep the plots aligned horizontally.
2. If you want to maintain the different y-axis scales for each plot, consider creating separate figures for each dataset and then combining them into a single figure using gridspec.
3. You can also use multiple y-axes with the same range but different scales for both plots.

Here's the updated code:

import matplotlib.pyplot as plt
import numpy as np

data1 = np.array([
 [-2, 22.8],
 [-2, 19.6],
 [-2, 0.3],
 [-1, 8.9],
 [-1, 13.7],
 [-1, 14.7],
 [ 0, 1.9],
 [ 0, -1.8],
 [ 1, -3],
 [ 1, -5.9],
 [ 1, -13.4],
 [ 2, -5.7],
 [ 2, -6.8],
]) 

data2 = np.array([
 [-2, 22.8],
 [-2, 0.3],
 [-2, 3.1],
 [-1, -1.7],
 [-1, 4.8],
 [-1, -0.7],
 [ 0, -2.6],
 [ 0, -0.03],
 [ 1, -5.7],
 [ 1, -1.5],
 [ 1, -3.9],
 [ 2, -21.5],
 [ 2, -7.7],
]) 


custom_annotations = ["K464E", "K472E", "R470E", "K464A", "M155E", "K472A", "M155A", "Q539A", "M155R", "D244A", "E247A", "E247R", "D244K"]


fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 5))

ax1.scatter(data1[:,0], data1[:,1])
ax1.set_title('Experimental Apo mHCN2 Channel')
ax1.set_xlabel(r'$\Delta q$')
ax1.set_ylabel(r'$\delta V_{0.5}$ Apo wild-type mHCN2 (mV)')

ax2.scatter(data2[:,0], data2[:,1])
ax2.set_title('Electrostatic Potential Calculations - Holo mHCN2 Channel')
ax2.set_xlabel(r'$\Delta q$')
ax2.set_ylabel(r'$\delta V_{0.5}$ Apo wild-type mHCN2 (mV)')

plt.axvline(0, c=(.5, .5, .5), ls= '--')
plt.axhline(0, c=(.5, .5, .5), ls= '--')



for i, txt in enumerate(custom_annotations):
 ax1.annotate(txt, (data1[i,0], data1[i,1]))
 ax2.annotate(txt, (data2[i,0], data2[i,1]))

ax1.set_ylim([-20, 25])
ax2.set_ylim([-20, 25])

plt.show()

I hope this helps!