I would use np.digitize
to do the bin sorting for you. This way you can easily apply any function and set the range you are interested in.
import numpy as np
import pylab as plt
N = 2000
total_bins = 10
# Sample data
X = np.random.random(size=N)*10
Y = X**2 + np.random.random(size=N)*X*10
bins = np.linspace(X.min(),X.max(), total_bins)
delta = bins[1]-bins[0]
idx = np.digitize(X,bins)
running_median = [np.median(Y[idx==k]) for k in range(total_bins)]
plt.scatter(X,Y,color='k',alpha=.2,s=2)
plt.plot(bins-delta/2,running_median,'r--',lw=4,alpha=.8)
plt.axis('tight')
plt.show()
As an example of the versatility of the method, let's add errorbars given by the standard deviation of each bin:
running_std = [Y[idx==k].std() for k in range(total_bins)]
plt.errorbar(bins-delta/2,running_median,
running_std,fmt=None)