spherical_kde
0.1.2
These packages allow you to do rudimentary kernel density estimation on a sphere. Suggestions for improvements/extensions welcome.
The fundamental principle is to replace the traditional Gaussian function used in kernel density estimation with the Von Mises-Fisher distribution.
Bandwidth estimation is still rough-and-ready.
import numpyfrom spherical_kde import SphericalKDEimport matplotlib.pyplot as pltimport cartopy.crsfrom matplotlib.gridspec import GridSpec, GridSpecFromSubplotSpec# Choose a seed for deterministic plotnumpy.random.seed(seed=0)# Set up a grid of figuresfig = plt.figure(figsize=(10, 10))gs_vert = GridSpec(3, 1)gs_lower = GridSpecFromSubplotSpec(1, 2, subplot_spec=gs_vert[1])fig.add_subplot(gs_vert[0], projection=cartopy.crs.Mollweide())fig.add_subplot(gs_lower[0], projection=cartopy.crs.Orthographic())fig.add_subplot(gs_lower[1], projection=cartopy.crs.Orthographic(-10, 45))fig.add_subplot(gs_vert[2], projection=cartopy.crs.PlateCarree())# Choose parameters for samplesnsamples = 100pi = numpy.pi# Generate some samples centered on (1,1) +/- 0.3 radianstheta_samples = numpy.random.normal(loc=1, scale=0.3, size=nsamples)phi_samples = numpy.random.normal(loc=1, scale=0.3, size=nsamples)phi_samples = numpy.mod(phi_samples, pi*2)kde_green = SphericalKDE(phi_samples, theta_samples)# Generate some samples centered on (-1,1) +/- 0.4 radianstheta_samples = numpy.random.normal(loc=1, scale=0.4, size=nsamples)phi_samples = numpy.random.normal(loc=-1, scale=0.4, size=nsamples)phi_samples = numpy.mod(phi_samples, pi*2)kde_red = SphericalKDE(phi_samples, theta_samples)# Generate a spread of samples along latitude 2, height 0.1theta_samples = numpy.random.normal(loc=2, scale=0.1, size=nsamples)phi_samples = numpy.random.uniform(low=-pi/2, high=pi/2, size=nsamples)phi_samples = numpy.mod(phi_samples, pi*2)kde_blue = SphericalKDE(phi_samples, theta_samples, bandwidth=0.1)for ax in fig.axes:ax.set_global()ax.gridlines()ax.coastlines(linewidth=0.1)kde_green.plot(ax, 'g')kde_green.plot_samples(ax)kde_red.plot(ax, 'r')kde_blue.plot(ax, 'b')# Save to plotfig.tight_layout()fig.savefig('plot.png') 
