I’m using the
pm.gp.Latent.prior for a 2D function on polar-like coordinates (rho, theta).
X support (shape (n_points, 2)) of coordinates features angles -pi < theta < pi as the second dimension. Therefore, I need a kernel in the theta dimension which finds the closest angle difference between two angles.
First I tried modifying the
ExpQuad kernel to always take the smallest angle like so
import pymc3 as pm import numpy as np import theano.tensor as tt _PI = np.pi _2PI = 2*_PI class StationaryAngleWrap(pm.gp.cov.Stationary): def square_dist(self, X, Xs): dist = tt.sqrt(super(StationaryAngleWrap, self).square_dist(X, Xs) + 1e-12) # calculate differences without squares return tt.square(tt.switch(dist > _PI, dist - _2PI, dist)) # get smallest possible angle class ExpQuadAngleWrap(StationaryAngleWrap, pm.gp.cov.ExpQuad): pass
However, this usually results in the model failing to initialize with bad initial energy. My guess is that
switch complicated the gradient too much.
Therefore, I tried creating a compound kernel which tries different angle difference offsets and should give high covariance for the smallest difference.
class StationaryOffset(pm.gp.cov.Stationary): def __init__(self, input_dim, ls=None, ls_inv=None, active_dims=None, offset=0): super(StationaryOffset, self).__init__(input_dim, ls, ls_inv, active_dims) self._offset = offset def square_dist(self, X, Xs): return super(StationaryOffset, self).square_dist(X - self._offset, Xs) class ExpQuadOffset(StationaryOffset, pm.gp.cov.ExpQuad): pass cov_theta = (ExpQuadOffset(2, theta_o_ls, active_dims=) + ExpQuadOffset(2, theta_o_ls, active_dims=, offset=2*np.pi) + ExpQuadOffset(2, theta_o_ls, active_dims=, offset=-2*np.pi))
This approach mostly works, but sometimes give strange results at the -pi ~ pi discontinuiti in the coordinates.
Any suggestions how to improve either of these solutions or a different solution is welcome.