mirror of
https://github.com/SheffieldML/GPy.git
synced 2026-06-02 14:45:15 +02:00
Merge branch 'master' of github.com:SheffieldML/GPy
Conflicts: GPy/examples/__init__.py
This commit is contained in:
commit
7c3c2fc9c0
46 changed files with 829 additions and 559 deletions
|
|
@ -9,3 +9,10 @@ import util
|
||||||
import examples
|
import examples
|
||||||
from core import priors
|
from core import priors
|
||||||
import likelihoods
|
import likelihoods
|
||||||
|
import testing
|
||||||
|
from numpy.testing import Tester
|
||||||
|
from nose.tools import nottest
|
||||||
|
|
||||||
|
@nottest
|
||||||
|
def tests():
|
||||||
|
Tester(testing).test(verbose=10)
|
||||||
|
|
|
||||||
|
|
@ -5,3 +5,4 @@ import classification
|
||||||
import regression
|
import regression
|
||||||
import dimensionality_reduction
|
import dimensionality_reduction
|
||||||
import non_gaussian
|
import non_gaussian
|
||||||
|
import tutorials
|
||||||
|
|
|
||||||
|
|
@ -1,56 +0,0 @@
|
||||||
# The detailed explanations of the commands used in this file can be found in the tutorial section
|
|
||||||
|
|
||||||
import pylab as pb
|
|
||||||
pb.ion()
|
|
||||||
import numpy as np
|
|
||||||
import GPy
|
|
||||||
|
|
||||||
X = np.random.uniform(-3.,3.,(20,1))
|
|
||||||
Y = np.sin(X) + np.random.randn(20,1)*0.05
|
|
||||||
|
|
||||||
kernel = GPy.kern.rbf(D=1, variance=1., lengthscale=1.)
|
|
||||||
|
|
||||||
m = GPy.models.GP_regression(X,Y,kernel)
|
|
||||||
|
|
||||||
print m
|
|
||||||
m.plot()
|
|
||||||
|
|
||||||
m.constrain_positive('')
|
|
||||||
|
|
||||||
m.unconstrain('') # Required to remove the previous constrains
|
|
||||||
m.constrain_positive('rbf_variance')
|
|
||||||
m.constrain_bounded('lengthscale',1.,10. )
|
|
||||||
m.constrain_fixed('noise',0.0025)
|
|
||||||
|
|
||||||
m.optimize()
|
|
||||||
|
|
||||||
m.optimize_restarts(Nrestarts = 10)
|
|
||||||
|
|
||||||
###########################
|
|
||||||
# 2-dimensional example #
|
|
||||||
###########################
|
|
||||||
|
|
||||||
import pylab as pb
|
|
||||||
pb.ion()
|
|
||||||
import numpy as np
|
|
||||||
import GPy
|
|
||||||
|
|
||||||
# sample inputs and outputs
|
|
||||||
X = np.random.uniform(-3.,3.,(50,2))
|
|
||||||
Y = np.sin(X[:,0:1]) * np.sin(X[:,1:2])+np.random.randn(50,1)*0.05
|
|
||||||
|
|
||||||
# define kernel
|
|
||||||
ker = GPy.kern.Matern52(2,ARD=True) + GPy.kern.white(2)
|
|
||||||
|
|
||||||
# create simple GP model
|
|
||||||
m = GPy.models.GP_regression(X,Y,ker)
|
|
||||||
|
|
||||||
# contrain all parameters to be positive
|
|
||||||
m.constrain_positive('')
|
|
||||||
|
|
||||||
# optimize and plot
|
|
||||||
pb.figure()
|
|
||||||
m.optimize('tnc', max_f_eval = 1000)
|
|
||||||
|
|
||||||
m.plot()
|
|
||||||
print(m)
|
|
||||||
|
|
@ -1,139 +0,0 @@
|
||||||
# The detailed explanations of the commands used in this file can be found in the tutorial section
|
|
||||||
|
|
||||||
import pylab as pb
|
|
||||||
import numpy as np
|
|
||||||
import GPy
|
|
||||||
pb.ion()
|
|
||||||
|
|
||||||
ker1 = GPy.kern.rbf(1) # Equivalent to ker1 = GPy.kern.rbf(D=1, variance=1., lengthscale=1.)
|
|
||||||
ker2 = GPy.kern.rbf(D=1, variance = .75, lengthscale=2.)
|
|
||||||
ker3 = GPy.kern.rbf(1, .5, .5)
|
|
||||||
|
|
||||||
print ker2
|
|
||||||
ker1.plot()
|
|
||||||
ker2.plot()
|
|
||||||
ker3.plot()
|
|
||||||
|
|
||||||
k1 = GPy.kern.rbf(1,1.,2.)
|
|
||||||
k2 = GPy.kern.Matern32(1, 0.5, 0.2)
|
|
||||||
|
|
||||||
# Product of kernels
|
|
||||||
k_prod = k1.prod(k2)
|
|
||||||
k_prodorth = k1.prod_orthogonal(k2)
|
|
||||||
|
|
||||||
# Sum of kernels
|
|
||||||
k_add = k1.add(k2)
|
|
||||||
k_addorth = k1.add_orthogonal(k2)
|
|
||||||
|
|
||||||
pb.figure(figsize=(8,8))
|
|
||||||
pb.subplot(2,2,1)
|
|
||||||
k_prod.plot()
|
|
||||||
pb.title('prod')
|
|
||||||
pb.subplot(2,2,2)
|
|
||||||
k_prodorth.plot()
|
|
||||||
pb.title('prod_orthogonal')
|
|
||||||
pb.subplot(2,2,3)
|
|
||||||
k_add.plot()
|
|
||||||
pb.title('add')
|
|
||||||
pb.subplot(2,2,4)
|
|
||||||
k_addorth.plot()
|
|
||||||
pb.title('add_orthogonal')
|
|
||||||
pb.subplots_adjust(wspace=0.3, hspace=0.3)
|
|
||||||
|
|
||||||
k1 = GPy.kern.rbf(1,1.,2)
|
|
||||||
k2 = GPy.kern.periodic_Matern52(1,variance=1e3, lengthscale=1, period = 1.5, lower=-5., upper = 5)
|
|
||||||
|
|
||||||
k = k1 * k2 # equivalent to k = k1.prod(k2)
|
|
||||||
print k
|
|
||||||
|
|
||||||
# Simulate sample paths
|
|
||||||
X = np.linspace(-5,5,501)[:,None]
|
|
||||||
Y = np.random.multivariate_normal(np.zeros(501),k.K(X),1)
|
|
||||||
|
|
||||||
# plot
|
|
||||||
pb.figure(figsize=(10,4))
|
|
||||||
pb.subplot(1,2,1)
|
|
||||||
k.plot()
|
|
||||||
pb.subplot(1,2,2)
|
|
||||||
pb.plot(X,Y.T)
|
|
||||||
pb.ylabel("Sample path")
|
|
||||||
pb.subplots_adjust(wspace=0.3)
|
|
||||||
|
|
||||||
k = (k1+k2)*(k1+k2)
|
|
||||||
print k.parts[0].name, '\n', k.parts[1].name, '\n', k.parts[2].name, '\n', k.parts[3].name
|
|
||||||
|
|
||||||
k1 = GPy.kern.rbf(1)
|
|
||||||
k2 = GPy.kern.Matern32(1)
|
|
||||||
k3 = GPy.kern.white(1)
|
|
||||||
|
|
||||||
k = k1 + k2 + k3
|
|
||||||
print k
|
|
||||||
|
|
||||||
k.constrain_positive('var')
|
|
||||||
k.constrain_fixed(np.array([1]),1.75)
|
|
||||||
k.tie_param('len')
|
|
||||||
k.unconstrain('white')
|
|
||||||
k.constrain_bounded('white',lower=1e-5,upper=.5)
|
|
||||||
print k
|
|
||||||
|
|
||||||
k_cst = GPy.kern.bias(1,variance=1.)
|
|
||||||
k_mat = GPy.kern.Matern52(1,variance=1., lengthscale=3)
|
|
||||||
Kanova = (k_cst + k_mat).prod_orthogonal(k_cst + k_mat)
|
|
||||||
print Kanova
|
|
||||||
|
|
||||||
# sample inputs and outputs
|
|
||||||
X = np.random.uniform(-3.,3.,(40,2))
|
|
||||||
Y = 0.5*X[:,:1] + 0.5*X[:,1:] + 2*np.sin(X[:,:1]) * np.sin(X[:,1:])
|
|
||||||
|
|
||||||
# Create GP regression model
|
|
||||||
m = GPy.models.GP_regression(X,Y,Kanova)
|
|
||||||
pb.figure(figsize=(5,5))
|
|
||||||
m.plot()
|
|
||||||
|
|
||||||
pb.figure(figsize=(20,3))
|
|
||||||
pb.subplots_adjust(wspace=0.5)
|
|
||||||
pb.subplot(1,5,1)
|
|
||||||
m.plot()
|
|
||||||
pb.subplot(1,5,2)
|
|
||||||
pb.ylabel("= ",rotation='horizontal',fontsize='30')
|
|
||||||
pb.subplot(1,5,3)
|
|
||||||
m.plot(which_functions=[False,True,False,False])
|
|
||||||
pb.ylabel("cst +",rotation='horizontal',fontsize='30')
|
|
||||||
pb.subplot(1,5,4)
|
|
||||||
m.plot(which_functions=[False,False,True,False])
|
|
||||||
pb.ylabel("+ ",rotation='horizontal',fontsize='30')
|
|
||||||
pb.subplot(1,5,5)
|
|
||||||
pb.ylabel("+ ",rotation='horizontal',fontsize='30')
|
|
||||||
m.plot(which_functions=[False,False,False,True])
|
|
||||||
|
|
||||||
import pylab as pb
|
|
||||||
import numpy as np
|
|
||||||
import GPy
|
|
||||||
pb.ion()
|
|
||||||
|
|
||||||
ker1 = GPy.kern.rbf(D=1) # Equivalent to ker1 = GPy.kern.rbf(D=1, variance=1., lengthscale=1.)
|
|
||||||
ker2 = GPy.kern.rbf(D=1, variance = .75, lengthscale=3.)
|
|
||||||
ker3 = GPy.kern.rbf(1, .5, .25)
|
|
||||||
|
|
||||||
ker1.plot()
|
|
||||||
ker2.plot()
|
|
||||||
ker3.plot()
|
|
||||||
#pb.savefig("Figures/tuto_kern_overview_basicdef.png")
|
|
||||||
|
|
||||||
kernels = [GPy.kern.rbf(1), GPy.kern.exponential(1), GPy.kern.Matern32(1), GPy.kern.Matern52(1), GPy.kern.Brownian(1), GPy.kern.bias(1), GPy.kern.linear(1), GPy.kern.spline(1), GPy.kern.periodic_exponential(1), GPy.kern.periodic_Matern32(1), GPy.kern.periodic_Matern52(1), GPy.kern.white(1)]
|
|
||||||
kernel_names = ["GPy.kern.rbf", "GPy.kern.exponential", "GPy.kern.Matern32", "GPy.kern.Matern52", "GPy.kern.Brownian", "GPy.kern.bias", "GPy.kern.linear", "GPy.kern.spline", "GPy.kern.periodic_exponential", "GPy.kern.periodic_Matern32", "GPy.kern.periodic_Matern52", "GPy.kern.white"]
|
|
||||||
|
|
||||||
pb.figure(figsize=(16,12))
|
|
||||||
pb.subplots_adjust(wspace=.5, hspace=.5)
|
|
||||||
for i, kern in enumerate(kernels):
|
|
||||||
pb.subplot(3,4,i+1)
|
|
||||||
kern.plot(x=7.5,plot_limits=[0.00001,15.])
|
|
||||||
pb.title(kernel_names[i]+ '\n')
|
|
||||||
|
|
||||||
# actual plot for the noise
|
|
||||||
i = 11
|
|
||||||
X = np.linspace(0.,15.,201)
|
|
||||||
WN = 0*X
|
|
||||||
WN[100] = 1.
|
|
||||||
pb.subplot(3,4,i+1)
|
|
||||||
pb.plot(X,WN,'b')
|
|
||||||
201
GPy/examples/tutorials.py
Normal file
201
GPy/examples/tutorials.py
Normal file
|
|
@ -0,0 +1,201 @@
|
||||||
|
# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
|
||||||
|
# Licensed under the BSD 3-clause license (see LICENSE.txt)
|
||||||
|
|
||||||
|
|
||||||
|
"""
|
||||||
|
Code of Tutorials
|
||||||
|
"""
|
||||||
|
|
||||||
|
def tuto_GP_regression():
|
||||||
|
"""The detailed explanations of the commands used in this file can be found in the tutorial section"""
|
||||||
|
|
||||||
|
import pylab as pb
|
||||||
|
pb.ion()
|
||||||
|
import numpy as np
|
||||||
|
import GPy
|
||||||
|
|
||||||
|
X = np.random.uniform(-3.,3.,(20,1))
|
||||||
|
Y = np.sin(X) + np.random.randn(20,1)*0.05
|
||||||
|
|
||||||
|
kernel = GPy.kern.rbf(D=1, variance=1., lengthscale=1.)
|
||||||
|
|
||||||
|
m = GPy.models.GP_regression(X,Y,kernel)
|
||||||
|
|
||||||
|
print m
|
||||||
|
m.plot()
|
||||||
|
|
||||||
|
m.constrain_positive('')
|
||||||
|
|
||||||
|
m.unconstrain('') # Required to remove the previous constrains
|
||||||
|
m.constrain_positive('rbf_variance')
|
||||||
|
m.constrain_bounded('lengthscale',1.,10. )
|
||||||
|
m.constrain_fixed('noise',0.0025)
|
||||||
|
|
||||||
|
m.optimize()
|
||||||
|
|
||||||
|
m.optimize_restarts(Nrestarts = 10)
|
||||||
|
|
||||||
|
###########################
|
||||||
|
# 2-dimensional example #
|
||||||
|
###########################
|
||||||
|
|
||||||
|
import pylab as pb
|
||||||
|
pb.ion()
|
||||||
|
import numpy as np
|
||||||
|
import GPy
|
||||||
|
|
||||||
|
# sample inputs and outputs
|
||||||
|
X = np.random.uniform(-3.,3.,(50,2))
|
||||||
|
Y = np.sin(X[:,0:1]) * np.sin(X[:,1:2])+np.random.randn(50,1)*0.05
|
||||||
|
|
||||||
|
# define kernel
|
||||||
|
ker = GPy.kern.Matern52(2,ARD=True) + GPy.kern.white(2)
|
||||||
|
|
||||||
|
# create simple GP model
|
||||||
|
m = GPy.models.GP_regression(X,Y,ker)
|
||||||
|
|
||||||
|
# contrain all parameters to be positive
|
||||||
|
m.constrain_positive('')
|
||||||
|
|
||||||
|
# optimize and plot
|
||||||
|
pb.figure()
|
||||||
|
m.optimize('tnc', max_f_eval = 1000)
|
||||||
|
|
||||||
|
m.plot()
|
||||||
|
print(m)
|
||||||
|
|
||||||
|
|
||||||
|
def tuto_kernel_overview():
|
||||||
|
"""The detailed explanations of the commands used in this file can be found in the tutorial section"""
|
||||||
|
import pylab as pb
|
||||||
|
import numpy as np
|
||||||
|
import GPy
|
||||||
|
pb.ion()
|
||||||
|
|
||||||
|
ker1 = GPy.kern.rbf(1) # Equivalent to ker1 = GPy.kern.rbf(D=1, variance=1., lengthscale=1.)
|
||||||
|
ker2 = GPy.kern.rbf(D=1, variance = .75, lengthscale=2.)
|
||||||
|
ker3 = GPy.kern.rbf(1, .5, .5)
|
||||||
|
|
||||||
|
print ker2
|
||||||
|
ker1.plot()
|
||||||
|
ker2.plot()
|
||||||
|
ker3.plot()
|
||||||
|
|
||||||
|
k1 = GPy.kern.rbf(1,1.,2.)
|
||||||
|
k2 = GPy.kern.Matern32(1, 0.5, 0.2)
|
||||||
|
|
||||||
|
# Product of kernels
|
||||||
|
k_prod = k1.prod(k2)
|
||||||
|
k_prodorth = k1.prod_orthogonal(k2)
|
||||||
|
|
||||||
|
# Sum of kernels
|
||||||
|
k_add = k1.add(k2)
|
||||||
|
k_addorth = k1.add_orthogonal(k2)
|
||||||
|
|
||||||
|
pb.figure(figsize=(8,8))
|
||||||
|
pb.subplot(2,2,1)
|
||||||
|
k_prod.plot()
|
||||||
|
pb.title('prod')
|
||||||
|
pb.subplot(2,2,2)
|
||||||
|
k_prodorth.plot()
|
||||||
|
pb.title('prod_orthogonal')
|
||||||
|
pb.subplot(2,2,3)
|
||||||
|
k_add.plot()
|
||||||
|
pb.title('add')
|
||||||
|
pb.subplot(2,2,4)
|
||||||
|
k_addorth.plot()
|
||||||
|
pb.title('add_orthogonal')
|
||||||
|
pb.subplots_adjust(wspace=0.3, hspace=0.3)
|
||||||
|
|
||||||
|
k1 = GPy.kern.rbf(1,1.,2)
|
||||||
|
k2 = GPy.kern.periodic_Matern52(1,variance=1e3, lengthscale=1, period = 1.5, lower=-5., upper = 5)
|
||||||
|
|
||||||
|
k = k1 * k2 # equivalent to k = k1.prod(k2)
|
||||||
|
print k
|
||||||
|
|
||||||
|
# Simulate sample paths
|
||||||
|
X = np.linspace(-5,5,501)[:,None]
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(501),k.K(X),1)
|
||||||
|
|
||||||
|
# plot
|
||||||
|
pb.figure(figsize=(10,4))
|
||||||
|
pb.subplot(1,2,1)
|
||||||
|
k.plot()
|
||||||
|
pb.subplot(1,2,2)
|
||||||
|
pb.plot(X,Y.T)
|
||||||
|
pb.ylabel("Sample path")
|
||||||
|
pb.subplots_adjust(wspace=0.3)
|
||||||
|
|
||||||
|
k = (k1+k2)*(k1+k2)
|
||||||
|
print k.parts[0].name, '\n', k.parts[1].name, '\n', k.parts[2].name, '\n', k.parts[3].name
|
||||||
|
|
||||||
|
k1 = GPy.kern.rbf(1)
|
||||||
|
k2 = GPy.kern.Matern32(1)
|
||||||
|
k3 = GPy.kern.white(1)
|
||||||
|
|
||||||
|
k = k1 + k2 + k3
|
||||||
|
print k
|
||||||
|
|
||||||
|
k.constrain_positive('var')
|
||||||
|
k.constrain_fixed(np.array([1]),1.75)
|
||||||
|
k.tie_param('len')
|
||||||
|
k.unconstrain('white')
|
||||||
|
k.constrain_bounded('white',lower=1e-5,upper=.5)
|
||||||
|
print k
|
||||||
|
|
||||||
|
k_cst = GPy.kern.bias(1,variance=1.)
|
||||||
|
k_mat = GPy.kern.Matern52(1,variance=1., lengthscale=3)
|
||||||
|
Kanova = (k_cst + k_mat).prod_orthogonal(k_cst + k_mat)
|
||||||
|
print Kanova
|
||||||
|
|
||||||
|
# sample inputs and outputs
|
||||||
|
X = np.random.uniform(-3.,3.,(40,2))
|
||||||
|
Y = 0.5*X[:,:1] + 0.5*X[:,1:] + 2*np.sin(X[:,:1]) * np.sin(X[:,1:])
|
||||||
|
|
||||||
|
# Create GP regression model
|
||||||
|
m = GPy.models.GP_regression(X,Y,Kanova)
|
||||||
|
pb.figure(figsize=(5,5))
|
||||||
|
m.plot()
|
||||||
|
|
||||||
|
pb.figure(figsize=(20,3))
|
||||||
|
pb.subplots_adjust(wspace=0.5)
|
||||||
|
pb.subplot(1,5,1)
|
||||||
|
m.plot()
|
||||||
|
pb.subplot(1,5,2)
|
||||||
|
pb.ylabel("= ",rotation='horizontal',fontsize='30')
|
||||||
|
pb.subplot(1,5,3)
|
||||||
|
m.plot(which_functions=[False,True,False,False])
|
||||||
|
pb.ylabel("cst +",rotation='horizontal',fontsize='30')
|
||||||
|
pb.subplot(1,5,4)
|
||||||
|
m.plot(which_functions=[False,False,True,False])
|
||||||
|
pb.ylabel("+ ",rotation='horizontal',fontsize='30')
|
||||||
|
pb.subplot(1,5,5)
|
||||||
|
pb.ylabel("+ ",rotation='horizontal',fontsize='30')
|
||||||
|
m.plot(which_functions=[False,False,False,True])
|
||||||
|
|
||||||
|
ker1 = GPy.kern.rbf(D=1) # Equivalent to ker1 = GPy.kern.rbf(D=1, variance=1., lengthscale=1.)
|
||||||
|
ker2 = GPy.kern.rbf(D=1, variance = .75, lengthscale=3.)
|
||||||
|
ker3 = GPy.kern.rbf(1, .5, .25)
|
||||||
|
|
||||||
|
ker1.plot()
|
||||||
|
ker2.plot()
|
||||||
|
ker3.plot()
|
||||||
|
#pb.savefig("Figures/tuto_kern_overview_basicdef.png")
|
||||||
|
|
||||||
|
kernels = [GPy.kern.rbf(1), GPy.kern.exponential(1), GPy.kern.Matern32(1), GPy.kern.Matern52(1), GPy.kern.Brownian(1), GPy.kern.bias(1), GPy.kern.linear(1), GPy.kern.spline(1), GPy.kern.periodic_exponential(1), GPy.kern.periodic_Matern32(1), GPy.kern.periodic_Matern52(1), GPy.kern.white(1)]
|
||||||
|
kernel_names = ["GPy.kern.rbf", "GPy.kern.exponential", "GPy.kern.Matern32", "GPy.kern.Matern52", "GPy.kern.Brownian", "GPy.kern.bias", "GPy.kern.linear", "GPy.kern.spline", "GPy.kern.periodic_exponential", "GPy.kern.periodic_Matern32", "GPy.kern.periodic_Matern52", "GPy.kern.white"]
|
||||||
|
|
||||||
|
pb.figure(figsize=(16,12))
|
||||||
|
pb.subplots_adjust(wspace=.5, hspace=.5)
|
||||||
|
for i, kern in enumerate(kernels):
|
||||||
|
pb.subplot(3,4,i+1)
|
||||||
|
kern.plot(x=7.5,plot_limits=[0.00001,15.])
|
||||||
|
pb.title(kernel_names[i]+ '\n')
|
||||||
|
|
||||||
|
# actual plot for the noise
|
||||||
|
i = 11
|
||||||
|
X = np.linspace(0.,15.,201)
|
||||||
|
WN = 0*X
|
||||||
|
WN[100] = 1.
|
||||||
|
pb.subplot(3,4,i+1)
|
||||||
|
pb.plot(X,WN,'b')
|
||||||
|
|
@ -76,7 +76,7 @@ class Matern32(kernpart):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
np.add(target,self.variance,target)
|
np.add(target,self.variance,target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters."""
|
"""derivative of the covariance matrix with respect to the parameters."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))
|
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))
|
||||||
|
|
@ -84,29 +84,29 @@ class Matern32(kernpart):
|
||||||
invdist = 1./np.where(dist!=0.,dist,np.inf)
|
invdist = 1./np.where(dist!=0.,dist,np.inf)
|
||||||
dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3
|
dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3
|
||||||
#dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
#dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
||||||
target[0] += np.sum(dvar*partial)
|
target[0] += np.sum(dvar*dL_dK)
|
||||||
if self.ARD == True:
|
if self.ARD == True:
|
||||||
dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
||||||
#dl = self.variance*dvar[:,:,None]*dist2M*invdist[:,:,None]
|
#dl = self.variance*dvar[:,:,None]*dist2M*invdist[:,:,None]
|
||||||
target[1:] += (dl*partial[:,:,None]).sum(0).sum(0)
|
target[1:] += (dl*dL_dK[:,:,None]).sum(0).sum(0)
|
||||||
else:
|
else:
|
||||||
dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist)) * dist2M.sum(-1)*invdist
|
dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist)) * dist2M.sum(-1)*invdist
|
||||||
#dl = self.variance*dvar*dist2M.sum(-1)*invdist
|
#dl = self.variance*dvar*dist2M.sum(-1)*invdist
|
||||||
target[1] += np.sum(dl*partial)
|
target[1] += np.sum(dl*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""derivative of the diagonal of the covariance matrix with respect to the parameters."""
|
"""derivative of the diagonal of the covariance matrix with respect to the parameters."""
|
||||||
target[0] += np.sum(partial)
|
target[0] += np.sum(dL_dKdiag)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to X."""
|
"""derivative of the covariance matrix with respect to X."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None]
|
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None]
|
||||||
ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf)
|
ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf)
|
||||||
dK_dX = - np.transpose(3*self.variance*dist*np.exp(-np.sqrt(3)*dist)*ddist_dX,(1,0,2))
|
dK_dX = - np.transpose(3*self.variance*dist*np.exp(-np.sqrt(3)*dist)*ddist_dX,(1,0,2))
|
||||||
target += np.sum(dK_dX*partial.T[:,:,None],0)
|
target += np.sum(dK_dX*dL_dK.T[:,:,None],0)
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def Gram_matrix(self,F,F1,F2,lower,upper):
|
def Gram_matrix(self,F,F1,F2,lower,upper):
|
||||||
|
|
|
||||||
|
|
@ -74,7 +74,7 @@ class Matern52(kernpart):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
np.add(target,self.variance,target)
|
np.add(target,self.variance,target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters."""
|
"""derivative of the covariance matrix with respect to the parameters."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))
|
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))
|
||||||
|
|
@ -82,29 +82,29 @@ class Matern52(kernpart):
|
||||||
dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3
|
dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3
|
||||||
dvar = (1+np.sqrt(5.)*dist+5./3*dist**2)*np.exp(-np.sqrt(5.)*dist)
|
dvar = (1+np.sqrt(5.)*dist+5./3*dist**2)*np.exp(-np.sqrt(5.)*dist)
|
||||||
dl = (self.variance * 5./3 * dist * (1 + np.sqrt(5.)*dist ) * np.exp(-np.sqrt(5.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
dl = (self.variance * 5./3 * dist * (1 + np.sqrt(5.)*dist ) * np.exp(-np.sqrt(5.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
||||||
target[0] += np.sum(dvar*partial)
|
target[0] += np.sum(dvar*dL_dK)
|
||||||
if self.ARD:
|
if self.ARD:
|
||||||
dl = (self.variance * 5./3 * dist * (1 + np.sqrt(5.)*dist ) * np.exp(-np.sqrt(5.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
dl = (self.variance * 5./3 * dist * (1 + np.sqrt(5.)*dist ) * np.exp(-np.sqrt(5.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
||||||
#dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
#dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist))[:,:,np.newaxis] * dist2M*invdist[:,:,np.newaxis]
|
||||||
target[1:] += (dl*partial[:,:,None]).sum(0).sum(0)
|
target[1:] += (dl*dL_dK[:,:,None]).sum(0).sum(0)
|
||||||
else:
|
else:
|
||||||
dl = (self.variance * 5./3 * dist * (1 + np.sqrt(5.)*dist ) * np.exp(-np.sqrt(5.)*dist)) * dist2M.sum(-1)*invdist
|
dl = (self.variance * 5./3 * dist * (1 + np.sqrt(5.)*dist ) * np.exp(-np.sqrt(5.)*dist)) * dist2M.sum(-1)*invdist
|
||||||
#dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist)) * dist2M.sum(-1)*invdist
|
#dl = (self.variance* 3 * dist * np.exp(-np.sqrt(3.)*dist)) * dist2M.sum(-1)*invdist
|
||||||
target[1] += np.sum(dl*partial)
|
target[1] += np.sum(dl*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""derivative of the diagonal of the covariance matrix with respect to the parameters."""
|
"""derivative of the diagonal of the covariance matrix with respect to the parameters."""
|
||||||
target[0] += np.sum(partial)
|
target[0] += np.sum(dL_dKdiag)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to X."""
|
"""derivative of the covariance matrix with respect to X."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None]
|
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None]
|
||||||
ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf)
|
ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf)
|
||||||
dK_dX = - np.transpose(self.variance*5./3*dist*(1+np.sqrt(5)*dist)*np.exp(-np.sqrt(5)*dist)*ddist_dX,(1,0,2))
|
dK_dX = - np.transpose(self.variance*5./3*dist*(1+np.sqrt(5)*dist)*np.exp(-np.sqrt(5)*dist)*ddist_dX,(1,0,2))
|
||||||
target += np.sum(dK_dX*partial.T[:,:,None],0)
|
target += np.sum(dK_dX*dL_dK.T[:,:,None],0)
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def Gram_matrix(self,F,F1,F2,F3,lower,upper):
|
def Gram_matrix(self,F,F1,F2,F3,lower,upper):
|
||||||
|
|
|
||||||
|
|
@ -35,16 +35,17 @@ class bias(kernpart):
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
target += self.variance
|
target += self.variance
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dKdiag,X,X2,target):
|
||||||
target += partial.sum()
|
target += dL_dKdiag.sum()
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
|
||||||
target += partial.sum()
|
|
||||||
|
|
||||||
def dK_dX(self, partial,X, X2, target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
|
target += dL_dKdiag.sum()
|
||||||
|
|
||||||
|
def dK_dX(self, dL_dK,X, X2, target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
#---------------------------------------#
|
#---------------------------------------#
|
||||||
|
|
@ -60,30 +61,29 @@ class bias(kernpart):
|
||||||
def psi2(self, Z, mu, S, target):
|
def psi2(self, Z, mu, S, target):
|
||||||
target += self.variance**2
|
target += self.variance**2
|
||||||
|
|
||||||
def dpsi0_dtheta(self, partial, Z, mu, S, target):
|
def dpsi0_dtheta(self, dL_dpsi0, Z, mu, S, target):
|
||||||
target += partial.sum()
|
target += dL_dpsi0.sum()
|
||||||
|
|
||||||
def dpsi1_dtheta(self, partial, Z, mu, S, target):
|
def dpsi1_dtheta(self, dL_dpsi1, Z, mu, S, target):
|
||||||
target += partial.sum()
|
target += dL_dpsi1.sum()
|
||||||
|
|
||||||
def dpsi2_dtheta(self, partial, Z, mu, S, target):
|
def dpsi2_dtheta(self, dL_dpsi2, Z, mu, S, target):
|
||||||
target += 2.*self.variance*partial.sum()
|
target += 2.*self.variance*dL_dpsi2.sum()
|
||||||
|
|
||||||
|
def dpsi0_dZ(self, dL_dpsi0, Z, mu, S, target):
|
||||||
def dpsi0_dZ(self, partial, Z, mu, S, target):
|
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi0_dmuS(self, partial, Z, mu, S, target_mu, target_S):
|
def dpsi0_dmuS(self, dL_dpsi0, Z, mu, S, target_mu, target_S):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi1_dZ(self, partial, Z, mu, S, target):
|
def dpsi1_dZ(self, dL_dpsi1, Z, mu, S, target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi1_dmuS(self, partial, Z, mu, S, target_mu, target_S):
|
def dpsi1_dmuS(self, dL_dpsi1, Z, mu, S, target_mu, target_S):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi2_dZ(self, partial, Z, mu, S, target):
|
def dpsi2_dZ(self, dL_dpsi2, Z, mu, S, target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi2_dmuS(self, partial, Z, mu, S, target_mu, target_S):
|
def dpsi2_dmuS(self, dL_dpsi2, Z, mu, S, target_mu, target_S):
|
||||||
pass
|
pass
|
||||||
|
|
|
||||||
|
|
@ -53,7 +53,7 @@ class coregionalise(kernpart):
|
||||||
def Kdiag(self,index,target):
|
def Kdiag(self,index,target):
|
||||||
target += np.diag(self.B)[np.asarray(index,dtype=np.int).flatten()]
|
target += np.diag(self.B)[np.asarray(index,dtype=np.int).flatten()]
|
||||||
|
|
||||||
def dK_dtheta(self,partial,index,index2,target):
|
def dK_dtheta(self,dL_dK,index,index2,target):
|
||||||
index = np.asarray(index,dtype=np.int)
|
index = np.asarray(index,dtype=np.int)
|
||||||
if index2 is None:
|
if index2 is None:
|
||||||
index2 = index
|
index2 = index
|
||||||
|
|
@ -62,28 +62,28 @@ class coregionalise(kernpart):
|
||||||
ii,jj = np.meshgrid(index,index2)
|
ii,jj = np.meshgrid(index,index2)
|
||||||
ii,jj = ii.T, jj.T
|
ii,jj = ii.T, jj.T
|
||||||
|
|
||||||
partial_small = np.zeros_like(self.B)
|
dL_dK_small = np.zeros_like(self.B)
|
||||||
for i in range(self.Nout):
|
for i in range(self.Nout):
|
||||||
for j in range(self.Nout):
|
for j in range(self.Nout):
|
||||||
tmp = np.sum(partial[(ii==i)*(jj==j)])
|
tmp = np.sum(dL_dK[(ii==i)*(jj==j)])
|
||||||
partial_small[i,j] = tmp
|
dL_dK_small[i,j] = tmp
|
||||||
|
|
||||||
dkappa = np.diag(partial_small)
|
dkappa = np.diag(dL_dK_small)
|
||||||
partial_small += partial_small.T
|
dL_dK_small += dL_dK_small.T
|
||||||
dW = (self.W[:,None,:]*partial_small[:,:,None]).sum(0)
|
dW = (self.W[:,None,:]*dL_dK_small[:,:,None]).sum(0)
|
||||||
|
|
||||||
target += np.hstack([dW.flatten(),dkappa])
|
target += np.hstack([dW.flatten(),dkappa])
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,index,target):
|
def dKdiag_dtheta(self,dL_dKdiag,index,target):
|
||||||
index = np.asarray(index,dtype=np.int).flatten()
|
index = np.asarray(index,dtype=np.int).flatten()
|
||||||
partial_small = np.zeros(self.Nout)
|
dL_dKdiag_small = np.zeros(self.Nout)
|
||||||
for i in range(self.Nout):
|
for i in range(self.Nout):
|
||||||
partial_small[i] += np.sum(partial[index==i])
|
dL_dKdiag_small[i] += np.sum(dL_dKdiag[index==i])
|
||||||
dW = 2.*self.W*partial_small[:,None]
|
dW = 2.*self.W*dL_dKdiag_small[:,None]
|
||||||
dkappa = partial_small
|
dkappa = dL_dKdiag_small
|
||||||
target += np.hstack([dW.flatten(),dkappa])
|
target += np.hstack([dW.flatten(),dkappa])
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -74,35 +74,35 @@ class exponential(kernpart):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
np.add(target,self.variance,target)
|
np.add(target,self.variance,target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters."""
|
"""derivative of the covariance matrix with respect to the parameters."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))
|
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))
|
||||||
invdist = 1./np.where(dist!=0.,dist,np.inf)
|
invdist = 1./np.where(dist!=0.,dist,np.inf)
|
||||||
dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3
|
dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3
|
||||||
dvar = np.exp(-dist)
|
dvar = np.exp(-dist)
|
||||||
target[0] += np.sum(dvar*partial)
|
target[0] += np.sum(dvar*dL_dK)
|
||||||
if self.ARD == True:
|
if self.ARD == True:
|
||||||
dl = self.variance*dvar[:,:,None]*dist2M*invdist[:,:,None]
|
dl = self.variance*dvar[:,:,None]*dist2M*invdist[:,:,None]
|
||||||
target[1:] += (dl*partial[:,:,None]).sum(0).sum(0)
|
target[1:] += (dl*dL_dK[:,:,None]).sum(0).sum(0)
|
||||||
else:
|
else:
|
||||||
dl = self.variance*dvar*dist2M.sum(-1)*invdist
|
dl = self.variance*dvar*dist2M.sum(-1)*invdist
|
||||||
target[1] += np.sum(dl*partial)
|
target[1] += np.sum(dl*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""derivative of the diagonal of the covariance matrix with respect to the parameters."""
|
"""derivative of the diagonal of the covariance matrix with respect to the parameters."""
|
||||||
#NB: derivative of diagonal elements wrt lengthscale is 0
|
#NB: derivative of diagonal elements wrt lengthscale is 0
|
||||||
target[0] += np.sum(partial)
|
target[0] += np.sum(dL_dKdiag)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to X."""
|
"""derivative of the covariance matrix with respect to X."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None]
|
dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None]
|
||||||
ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf)
|
ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf)
|
||||||
dK_dX = - np.transpose(self.variance*np.exp(-dist)*ddist_dX,(1,0,2))
|
dK_dX = - np.transpose(self.variance*np.exp(-dist)*ddist_dX,(1,0,2))
|
||||||
target += np.sum(dK_dX*partial.T[:,:,None],0)
|
target += np.sum(dK_dX*dL_dK.T[:,:,None],0)
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def Gram_matrix(self,F,F1,lower,upper):
|
def Gram_matrix(self,F,F1,lower,upper):
|
||||||
|
|
|
||||||
|
|
@ -271,10 +271,10 @@ class kern(parameterised):
|
||||||
[p.K(X[s1,i_s],X2[s2,i_s],target=target[s1,s2]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
[p.K(X[s1,i_s],X2[s2,i_s],target=target[s1,s2]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2=None,slices1=None,slices2=None):
|
def dK_dtheta(self,dL_dK,X,X2=None,slices1=None,slices2=None):
|
||||||
"""
|
"""
|
||||||
:param partial: An array of partial derivaties, dL_dK
|
:param dL_dK: An array of dL_dK derivaties, dL_dK
|
||||||
:type partial: Np.ndarray (N x M)
|
:type dL_dK: Np.ndarray (N x M)
|
||||||
:param X: Observed data inputs
|
:param X: Observed data inputs
|
||||||
:type X: np.ndarray (N x D)
|
:type X: np.ndarray (N x D)
|
||||||
:param X2: Observed dara inputs (optional, defaults to X)
|
:param X2: Observed dara inputs (optional, defaults to X)
|
||||||
|
|
@ -288,16 +288,16 @@ class kern(parameterised):
|
||||||
if X2 is None:
|
if X2 is None:
|
||||||
X2 = X
|
X2 = X
|
||||||
target = np.zeros(self.Nparam)
|
target = np.zeros(self.Nparam)
|
||||||
[p.dK_dtheta(partial[s1,s2],X[s1,i_s],X2[s2,i_s],target[ps]) for p,i_s,ps,s1,s2 in zip(self.parts, self.input_slices, self.param_slices, slices1, slices2)]
|
[p.dK_dtheta(dL_dK[s1,s2],X[s1,i_s],X2[s2,i_s],target[ps]) for p,i_s,ps,s1,s2 in zip(self.parts, self.input_slices, self.param_slices, slices1, slices2)]
|
||||||
|
|
||||||
return self._transform_gradients(target)
|
return self._transform_gradients(target)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2=None,slices1=None,slices2=None):
|
def dK_dX(self,dL_dK,X,X2=None,slices1=None,slices2=None):
|
||||||
if X2 is None:
|
if X2 is None:
|
||||||
X2 = X
|
X2 = X
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target = np.zeros_like(X)
|
target = np.zeros_like(X)
|
||||||
[p.dK_dX(partial[s1,s2],X[s1,i_s],X2[s2,i_s],target[s1,i_s]) for p, i_s, s1, s2 in zip(self.parts, self.input_slices, slices1, slices2)]
|
[p.dK_dX(dL_dK[s1,s2],X[s1,i_s],X2[s2,i_s],target[s1,i_s]) for p, i_s, s1, s2 in zip(self.parts, self.input_slices, slices1, slices2)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def Kdiag(self,X,slices=None):
|
def Kdiag(self,X,slices=None):
|
||||||
|
|
@ -307,20 +307,20 @@ class kern(parameterised):
|
||||||
[p.Kdiag(X[s,i_s],target=target[s]) for p,i_s,s in zip(self.parts,self.input_slices,slices)]
|
[p.Kdiag(X[s,i_s],target=target[s]) for p,i_s,s in zip(self.parts,self.input_slices,slices)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,slices=None):
|
def dKdiag_dtheta(self,dL_dKdiag,X,slices=None):
|
||||||
assert X.shape[1]==self.D
|
assert X.shape[1]==self.D
|
||||||
assert len(partial.shape)==1
|
assert len(dL_dKdiag.shape)==1
|
||||||
assert partial.size==X.shape[0]
|
assert dL_dKdiag.size==X.shape[0]
|
||||||
slices = self._process_slices(slices,False)
|
slices = self._process_slices(slices,False)
|
||||||
target = np.zeros(self.Nparam)
|
target = np.zeros(self.Nparam)
|
||||||
[p.dKdiag_dtheta(partial[s],X[s,i_s],target[ps]) for p,i_s,s,ps in zip(self.parts,self.input_slices,slices,self.param_slices)]
|
[p.dKdiag_dtheta(dL_dKdiag[s],X[s,i_s],target[ps]) for p,i_s,s,ps in zip(self.parts,self.input_slices,slices,self.param_slices)]
|
||||||
return self._transform_gradients(target)
|
return self._transform_gradients(target)
|
||||||
|
|
||||||
def dKdiag_dX(self, partial, X, slices=None):
|
def dKdiag_dX(self, dL_dKdiag, X, slices=None):
|
||||||
assert X.shape[1]==self.D
|
assert X.shape[1]==self.D
|
||||||
slices = self._process_slices(slices,False)
|
slices = self._process_slices(slices,False)
|
||||||
target = np.zeros_like(X)
|
target = np.zeros_like(X)
|
||||||
[p.dKdiag_dX(partial[s],X[s,i_s],target[s,i_s]) for p,i_s,s in zip(self.parts,self.input_slices,slices)]
|
[p.dKdiag_dX(dL_dKdiag[s],X[s,i_s],target[s,i_s]) for p,i_s,s in zip(self.parts,self.input_slices,slices)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def psi0(self,Z,mu,S,slices=None):
|
def psi0(self,Z,mu,S,slices=None):
|
||||||
|
|
@ -329,16 +329,16 @@ class kern(parameterised):
|
||||||
[p.psi0(Z,mu[s],S[s],target[s]) for p,s in zip(self.parts,slices)]
|
[p.psi0(Z,mu[s],S[s],target[s]) for p,s in zip(self.parts,slices)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dpsi0_dtheta(self,partial,Z,mu,S,slices=None):
|
def dpsi0_dtheta(self,dL_dpsi0,Z,mu,S,slices=None):
|
||||||
slices = self._process_slices(slices,False)
|
slices = self._process_slices(slices,False)
|
||||||
target = np.zeros(self.Nparam)
|
target = np.zeros(self.Nparam)
|
||||||
[p.dpsi0_dtheta(partial[s],Z,mu[s],S[s],target[ps]) for p,ps,s in zip(self.parts, self.param_slices,slices)]
|
[p.dpsi0_dtheta(dL_dpsi0[s],Z,mu[s],S[s],target[ps]) for p,ps,s in zip(self.parts, self.param_slices,slices)]
|
||||||
return self._transform_gradients(target)
|
return self._transform_gradients(target)
|
||||||
|
|
||||||
def dpsi0_dmuS(self,partial,Z,mu,S,slices=None):
|
def dpsi0_dmuS(self,dL_dpsi0,Z,mu,S,slices=None):
|
||||||
slices = self._process_slices(slices,False)
|
slices = self._process_slices(slices,False)
|
||||||
target_mu,target_S = np.zeros_like(mu),np.zeros_like(S)
|
target_mu,target_S = np.zeros_like(mu),np.zeros_like(S)
|
||||||
[p.dpsi0_dmuS(partial,Z,mu[s],S[s],target_mu[s],target_S[s]) for p,s in zip(self.parts,slices)]
|
[p.dpsi0_dmuS(dL_dpsi0,Z,mu[s],S[s],target_mu[s],target_S[s]) for p,s in zip(self.parts,slices)]
|
||||||
return target_mu,target_S
|
return target_mu,target_S
|
||||||
|
|
||||||
def psi1(self,Z,mu,S,slices1=None,slices2=None):
|
def psi1(self,Z,mu,S,slices1=None,slices2=None):
|
||||||
|
|
@ -348,25 +348,25 @@ class kern(parameterised):
|
||||||
[p.psi1(Z[s2],mu[s1],S[s1],target[s1,s2]) for p,s1,s2 in zip(self.parts,slices1,slices2)]
|
[p.psi1(Z[s2],mu[s1],S[s1],target[s1,s2]) for p,s1,s2 in zip(self.parts,slices1,slices2)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dpsi1_dtheta(self,partial,Z,mu,S,slices1=None,slices2=None):
|
def dpsi1_dtheta(self,dL_dpsi1,Z,mu,S,slices1=None,slices2=None):
|
||||||
"""N,M,(Ntheta)"""
|
"""N,M,(Ntheta)"""
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target = np.zeros((self.Nparam))
|
target = np.zeros((self.Nparam))
|
||||||
[p.dpsi1_dtheta(partial[s2,s1],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[ps]) for p,ps,s1,s2,i_s in zip(self.parts, self.param_slices,slices1,slices2,self.input_slices)]
|
[p.dpsi1_dtheta(dL_dpsi1[s2,s1],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[ps]) for p,ps,s1,s2,i_s in zip(self.parts, self.param_slices,slices1,slices2,self.input_slices)]
|
||||||
return self._transform_gradients(target)
|
return self._transform_gradients(target)
|
||||||
|
|
||||||
def dpsi1_dZ(self,partial,Z,mu,S,slices1=None,slices2=None):
|
def dpsi1_dZ(self,dL_dpsi1,Z,mu,S,slices1=None,slices2=None):
|
||||||
"""N,M,Q"""
|
"""N,M,Q"""
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target = np.zeros_like(Z)
|
target = np.zeros_like(Z)
|
||||||
[p.dpsi1_dZ(partial[s2,s1],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[s2,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
[p.dpsi1_dZ(dL_dpsi1[s2,s1],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[s2,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dpsi1_dmuS(self,partial,Z,mu,S,slices1=None,slices2=None):
|
def dpsi1_dmuS(self,dL_dpsi1,Z,mu,S,slices1=None,slices2=None):
|
||||||
"""return shapes are N,M,Q"""
|
"""return shapes are N,M,Q"""
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target_mu, target_S = np.zeros((2,mu.shape[0],mu.shape[1]))
|
target_mu, target_S = np.zeros((2,mu.shape[0],mu.shape[1]))
|
||||||
[p.dpsi1_dmuS(partial[s2,s1],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target_mu[s1,i_s],target_S[s1,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
[p.dpsi1_dmuS(dL_dpsi1[s2,s1],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target_mu[s1,i_s],target_S[s1,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
||||||
return target_mu, target_S
|
return target_mu, target_S
|
||||||
|
|
||||||
def psi2(self,Z,mu,S,slices1=None,slices2=None):
|
def psi2(self,Z,mu,S,slices1=None,slices2=None):
|
||||||
|
|
@ -399,10 +399,11 @@ class kern(parameterised):
|
||||||
|
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dpsi2_dtheta(self,partial,partial1,Z,mu,S,slices1=None,slices2=None):
|
def dpsi2_dtheta(self,dL_dpsi2,Z,mu,S,slices1=None,slices2=None):
|
||||||
|
"""Returns shape (N,M,M,Ntheta)"""
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target = np.zeros(self.Nparam)
|
target = np.zeros(self.Nparam)
|
||||||
[p.dpsi2_dtheta(partial[s1,s2,s2],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[ps]) for p,i_s,s1,s2,ps in zip(self.parts,self.input_slices,slices1,slices2,self.param_slices)]
|
[p.dpsi2_dtheta(dL_dpsi2[s1,s2,s2],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[ps]) for p,i_s,s1,s2,ps in zip(self.parts,self.input_slices,slices1,slices2,self.param_slices)]
|
||||||
|
|
||||||
#compute the "cross" terms
|
#compute the "cross" terms
|
||||||
#TODO: better looping
|
#TODO: better looping
|
||||||
|
|
@ -416,11 +417,11 @@ class kern(parameterised):
|
||||||
pass
|
pass
|
||||||
#rbf X bias
|
#rbf X bias
|
||||||
elif p1.name=='bias' and p2.name=='rbf':
|
elif p1.name=='bias' and p2.name=='rbf':
|
||||||
p2.dpsi1_dtheta(partial.sum(1)*p1.variance,Z,mu,S,target[ps2])
|
p2.dpsi1_dtheta(dL_dpsi2.sum(1)*p1.variance,Z,mu,S,target[ps2])
|
||||||
p1.dpsi1_dtheta(partial.sum(1)*p2._psi1,Z,mu,S,target[ps1])
|
p1.dpsi1_dtheta(dL_dpsi2.sum(1)*p2._psi1,Z,mu,S,target[ps1])
|
||||||
elif p2.name=='bias' and p1.name=='rbf':
|
elif p2.name=='bias' and p1.name=='rbf':
|
||||||
p1.dpsi1_dtheta(partial.sum(1)*p2.variance,Z,mu,S,target[ps1])
|
p1.dpsi1_dtheta(dL_dpsi2.sum(1)*p2.variance,Z,mu,S,target[ps1])
|
||||||
p2.dpsi1_dtheta(partial.sum(1)*p1._psi1,Z,mu,S,target[ps2])
|
p2.dpsi1_dtheta(dL_dpsi2.sum(1)*p1._psi1,Z,mu,S,target[ps2])
|
||||||
#rbf X linear
|
#rbf X linear
|
||||||
elif p1.name=='linear' and p2.name=='rbf':
|
elif p1.name=='linear' and p2.name=='rbf':
|
||||||
raise NotImplementedError #TODO
|
raise NotImplementedError #TODO
|
||||||
|
|
@ -431,10 +432,10 @@ class kern(parameterised):
|
||||||
|
|
||||||
return self._transform_gradients(target)
|
return self._transform_gradients(target)
|
||||||
|
|
||||||
def dpsi2_dZ(self,partial,Z,mu,S,slices1=None,slices2=None):
|
def dpsi2_dZ(self,dL_dpsi2,Z,mu,S,slices1=None,slices2=None):
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target = np.zeros_like(Z)
|
target = np.zeros_like(Z)
|
||||||
[p.dpsi2_dZ(partial[s1,s2,s2],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[s2,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
[p.dpsi2_dZ(dL_dpsi2[s1,s2,s2],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target[s2,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
||||||
|
|
||||||
#compute the "cross" terms
|
#compute the "cross" terms
|
||||||
for p1, p2 in itertools.combinations(self.parts,2):
|
for p1, p2 in itertools.combinations(self.parts,2):
|
||||||
|
|
@ -443,9 +444,9 @@ class kern(parameterised):
|
||||||
pass
|
pass
|
||||||
#rbf X bias
|
#rbf X bias
|
||||||
elif p1.name=='bias' and p2.name=='rbf':
|
elif p1.name=='bias' and p2.name=='rbf':
|
||||||
target += p2.dpsi1_dX(partial.sum(1)*p1.variance,Z,mu,S,target)
|
p2.dpsi1_dX(dL_dpsi2.sum(1)*p1.variance,Z,mu,S,target)
|
||||||
elif p2.name=='bias' and p1.name=='rbf':
|
elif p2.name=='bias' and p1.name=='rbf':
|
||||||
target += p1.dpsi1_dZ(partial.sum(2)*p2.variance,Z,mu,S,target)
|
p1.dpsi1_dZ(dL_dpsi2.sum(2)*p2.variance,Z,mu,S,target)
|
||||||
#rbf X linear
|
#rbf X linear
|
||||||
elif p1.name=='linear' and p2.name=='rbf':
|
elif p1.name=='linear' and p2.name=='rbf':
|
||||||
raise NotImplementedError #TODO
|
raise NotImplementedError #TODO
|
||||||
|
|
@ -457,11 +458,11 @@ class kern(parameterised):
|
||||||
|
|
||||||
return target
|
return target
|
||||||
|
|
||||||
def dpsi2_dmuS(self,partial,Z,mu,S,slices1=None,slices2=None):
|
def dpsi2_dmuS(self,dL_dpsi2,Z,mu,S,slices1=None,slices2=None):
|
||||||
"""return shapes are N,M,M,Q"""
|
"""return shapes are N,M,M,Q"""
|
||||||
slices1, slices2 = self._process_slices(slices1,slices2)
|
slices1, slices2 = self._process_slices(slices1,slices2)
|
||||||
target_mu, target_S = np.zeros((2,mu.shape[0],mu.shape[1]))
|
target_mu, target_S = np.zeros((2,mu.shape[0],mu.shape[1]))
|
||||||
[p.dpsi2_dmuS(partial[s1,s2,s2],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target_mu[s1,i_s],target_S[s1,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
[p.dpsi2_dmuS(dL_dpsi2[s1,s2,s2],Z[s2,i_s],mu[s1,i_s],S[s1,i_s],target_mu[s1,i_s],target_S[s1,i_s]) for p,i_s,s1,s2 in zip(self.parts,self.input_slices,slices1,slices2)]
|
||||||
|
|
||||||
#compute the "cross" terms
|
#compute the "cross" terms
|
||||||
for p1, p2 in itertools.combinations(self.parts,2):
|
for p1, p2 in itertools.combinations(self.parts,2):
|
||||||
|
|
@ -470,9 +471,9 @@ class kern(parameterised):
|
||||||
pass
|
pass
|
||||||
#rbf X bias
|
#rbf X bias
|
||||||
elif p1.name=='bias' and p2.name=='rbf':
|
elif p1.name=='bias' and p2.name=='rbf':
|
||||||
target += p2.dpsi1_dmuS(partial.sum(1)*p1.variance,Z,mu,S,target_mu,target_S)
|
p2.dpsi1_dmuS(partial.sum(1)*p1.variance,Z,mu,S,target_mu,target_S)
|
||||||
elif p2.name=='bias' and p1.name=='rbf':
|
elif p2.name=='bias' and p1.name=='rbf':
|
||||||
target += p1.dpsi1_dmuS(partial.sum(2)*p2.variance,Z,mu,S,target_mu,target_S)
|
p1.dpsi1_dmuS(partial.sum(2)*p2.variance,Z,mu,S,target_mu,target_S)
|
||||||
#rbf X linear
|
#rbf X linear
|
||||||
elif p1.name=='linear' and p2.name=='rbf':
|
elif p1.name=='linear' and p2.name=='rbf':
|
||||||
raise NotImplementedError #TODO
|
raise NotImplementedError #TODO
|
||||||
|
|
|
||||||
|
|
@ -26,31 +26,31 @@ class kernpart(object):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def psi0(self,Z,mu,S,target):
|
def psi0(self,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi0_dtheta(self,partial,Z,mu,S,target):
|
def dpsi0_dtheta(self,dL_dpsi0,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi0_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi0_dmuS(self,dL_dpsi0,Z,mu,S,target_mu,target_S):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def psi1(self,Z,mu,S,target):
|
def psi1(self,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi1_dtheta(self,Z,mu,S,target):
|
def dpsi1_dtheta(self,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi1_dZ(self,partial,Z,mu,S,target):
|
def dpsi1_dZ(self,dL_dpsi1,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi1_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi1_dmuS(self,dL_dpsi1,Z,mu,S,target_mu,target_S):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def psi2(self,Z,mu,S,target):
|
def psi2(self,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi2_dZ(self,partial,Z,mu,S,target):
|
def dpsi2_dZ(self,dL_dpsi2,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi2_dtheta(self,partial,Z,mu,S,target):
|
def dpsi2_dtheta(self,dL_dpsi2,Z,mu,S,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dpsi2_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi2_dmuS(self,dL_dpsi2,Z,mu,S,target_mu,target_S):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
def dK_dX(self,X,X2,target):
|
def dK_dX(self,X,X2,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
|
|
|
||||||
|
|
@ -73,16 +73,16 @@ class linear(kernpart):
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
np.add(target,np.sum(self.variances*np.square(X),-1),target)
|
np.add(target,np.sum(self.variances*np.square(X),-1),target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
if self.ARD:
|
if self.ARD:
|
||||||
product = X[:,None,:]*X2[None,:,:]
|
product = X[:,None,:]*X2[None,:,:]
|
||||||
target += (partial[:,:,None]*product).sum(0).sum(0)
|
target += (dL_dK[:,:,None]*product).sum(0).sum(0)
|
||||||
else:
|
else:
|
||||||
self._K_computations(X, X2)
|
self._K_computations(X, X2)
|
||||||
target += np.sum(self._dot_product*partial)
|
target += np.sum(self._dot_product*dL_dK)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
target += (((X2[:, None, :] * self.variances)) * partial[:,:, None]).sum(0)
|
target += (((X2[:, None, :] * self.variances)) * dL_dK[:,:, None]).sum(0)
|
||||||
|
|
||||||
#---------------------------------------#
|
#---------------------------------------#
|
||||||
# PSI statistics #
|
# PSI statistics #
|
||||||
|
|
@ -92,40 +92,40 @@ class linear(kernpart):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
target += np.sum(self.variances*self.mu2_S,1)
|
target += np.sum(self.variances*self.mu2_S,1)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial, X, target):
|
def dKdiag_dtheta(self,dL_dKdiag, X, target):
|
||||||
tmp = partial[:,None]*X**2
|
tmp = dL_dKdiag[:,None]*X**2
|
||||||
if self.ARD:
|
if self.ARD:
|
||||||
target += tmp.sum(0)
|
target += tmp.sum(0)
|
||||||
else:
|
else:
|
||||||
target += tmp.sum()
|
target += tmp.sum()
|
||||||
|
|
||||||
def dpsi0_dtheta(self,partial,Z,mu,S,target):
|
def dpsi0_dtheta(self,dL_dpsi0,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
tmp = partial[:, None] * self.mu2_S
|
tmp = dL_dpsi0[:, None] * self.mu2_S
|
||||||
if self.ARD:
|
if self.ARD:
|
||||||
target += tmp.sum(0)
|
target += tmp.sum(0)
|
||||||
else:
|
else:
|
||||||
target += tmp.sum()
|
target += tmp.sum()
|
||||||
|
|
||||||
def dpsi0_dmuS(self,partial, Z,mu,S,target_mu,target_S):
|
def dpsi0_dmuS(self,dL_dpsi0, Z,mu,S,target_mu,target_S):
|
||||||
target_mu += partial[:, None] * (2.0*mu*self.variances)
|
target_mu += dL_dpsi0[:, None] * (2.0*mu*self.variances)
|
||||||
target_S += partial[:, None] * self.variances
|
target_S += dL_dpsi0[:, None] * self.variances
|
||||||
|
|
||||||
def psi1(self,Z,mu,S,target):
|
def psi1(self,Z,mu,S,target):
|
||||||
"""the variance, it does nothing"""
|
"""the variance, it does nothing"""
|
||||||
self.K(mu,Z,target)
|
self.K(mu,Z,target)
|
||||||
|
|
||||||
def dpsi1_dtheta(self,partial,Z,mu,S,target):
|
def dpsi1_dtheta(self,dL_dpsi1,Z,mu,S,target):
|
||||||
"""the variance, it does nothing"""
|
"""the variance, it does nothing"""
|
||||||
self.dK_dtheta(partial,mu,Z,target)
|
self.dK_dtheta(dL_dpsi1,mu,Z,target)
|
||||||
|
|
||||||
def dpsi1_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi1_dmuS(self,dL_dpsi1,Z,mu,S,target_mu,target_S):
|
||||||
"""Do nothing for S, it does not affect psi1"""
|
"""Do nothing for S, it does not affect psi1"""
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
target_mu += (partial.T[:,:, None]*(Z*self.variances)).sum(1)
|
target_mu += (dL_dpsi1.T[:,:, None]*(Z*self.variances)).sum(1)
|
||||||
|
|
||||||
def dpsi1_dZ(self,partial,Z,mu,S,target):
|
def dpsi1_dZ(self,dL_dpsi1,Z,mu,S,target):
|
||||||
self.dK_dX(partial.T,Z,mu,target)
|
self.dK_dX(dL_dpsi1.T,Z,mu,target)
|
||||||
|
|
||||||
def psi2(self,Z,mu,S,target):
|
def psi2(self,Z,mu,S,target):
|
||||||
"""
|
"""
|
||||||
|
|
@ -135,25 +135,25 @@ class linear(kernpart):
|
||||||
psi2 = self.ZZ*np.square(self.variances)*self.mu2_S[:, None, None, :]
|
psi2 = self.ZZ*np.square(self.variances)*self.mu2_S[:, None, None, :]
|
||||||
target += psi2.sum(-1)
|
target += psi2.sum(-1)
|
||||||
|
|
||||||
def dpsi2_dtheta(self,partial,Z,mu,S,target):
|
def dpsi2_dtheta(self,dL_dpsi2,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
tmp = (partial[:,:,:,None]*(2.*self.ZZ*self.mu2_S[:,None,None,:]*self.variances))
|
tmp = (dL_dpsi2[:,:,:,None]*(2.*self.ZZ*self.mu2_S[:,None,None,:]*self.variances))
|
||||||
if self.ARD:
|
if self.ARD:
|
||||||
target += tmp.sum(0).sum(0).sum(0)
|
target += tmp.sum(0).sum(0).sum(0)
|
||||||
else:
|
else:
|
||||||
target += tmp.sum()
|
target += tmp.sum()
|
||||||
|
|
||||||
def dpsi2_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi2_dmuS(self,dL_dpsi2,Z,mu,S,target_mu,target_S):
|
||||||
"""Think N,M,M,Q """
|
"""Think N,M,M,Q """
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
tmp = self.ZZ*np.square(self.variances) # M,M,Q
|
tmp = self.ZZ*np.square(self.variances) # M,M,Q
|
||||||
target_mu += (partial[:,:,:,None]*tmp*2.*mu[:,None,None,:]).sum(1).sum(1)
|
target_mu += (dL_dpsi2[:,:,:,None]*tmp*2.*mu[:,None,None,:]).sum(1).sum(1)
|
||||||
target_S += (partial[:,:,:,None]*tmp).sum(1).sum(1)
|
target_S += (dL_dpsi2[:,:,:,None]*tmp).sum(1).sum(1)
|
||||||
|
|
||||||
def dpsi2_dZ(self,partial,Z,mu,S,target):
|
def dpsi2_dZ(self,dL_dpsi2,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
mu2_S = np.sum(self.mu2_S,0)# Q,
|
mu2_S = np.sum(self.mu2_S,0)# Q,
|
||||||
target += (partial[:,:,:,None] * (self.mu2_S[:,None,None,:]*(Z*np.square(self.variances)[None,:])[None,None,:,:])).sum(0).sum(1)
|
target += (dL_dpsi2[:,:,:,None] * (self.mu2_S[:,None,None,:]*(Z*np.square(self.variances)[None,:])[None,None,:,:])).sum(0).sum(1)
|
||||||
|
|
||||||
#---------------------------------------#
|
#---------------------------------------#
|
||||||
# Precomputations #
|
# Precomputations #
|
||||||
|
|
|
||||||
|
|
@ -101,7 +101,7 @@ class periodic_Matern32(kernpart):
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
np.add(target,np.diag(mdot(FX,self.Gi,FX.T)),target)
|
np.add(target,np.diag(mdot(FX,self.Gi,FX.T)),target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters (shape is NxMxNparam)"""
|
"""derivative of the covariance matrix with respect to the parameters (shape is NxMxNparam)"""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
|
|
@ -166,13 +166,13 @@ class periodic_Matern32(kernpart):
|
||||||
dK_dper = mdot(dFX_dper,self.Gi,FX2.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX2.T) + mdot(FX,self.Gi,dFX2_dper.T)
|
dK_dper = mdot(dFX_dper,self.Gi,FX2.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX2.T) + mdot(FX,self.Gi,dFX2_dper.T)
|
||||||
|
|
||||||
# np.add(target[:,:,0],dK_dvar, target[:,:,0])
|
# np.add(target[:,:,0],dK_dvar, target[:,:,0])
|
||||||
target[0] += np.sum(dK_dvar*partial)
|
target[0] += np.sum(dK_dvar*dL_dK)
|
||||||
#np.add(target[:,:,1],dK_dlen, target[:,:,1])
|
#np.add(target[:,:,1],dK_dlen, target[:,:,1])
|
||||||
target[1] += np.sum(dK_dlen*partial)
|
target[1] += np.sum(dK_dlen*dL_dK)
|
||||||
#np.add(target[:,:,2],dK_dper, target[:,:,2])
|
#np.add(target[:,:,2],dK_dper, target[:,:,2])
|
||||||
target[2] += np.sum(dK_dper*partial)
|
target[2] += np.sum(dK_dper*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""derivative of the diagonal covariance matrix with respect to the parameters"""
|
"""derivative of the diagonal covariance matrix with respect to the parameters"""
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
|
|
||||||
|
|
@ -231,6 +231,6 @@ class periodic_Matern32(kernpart):
|
||||||
|
|
||||||
dK_dper = 2* mdot(dFX_dper,self.Gi,FX.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX.T)
|
dK_dper = 2* mdot(dFX_dper,self.Gi,FX.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX.T)
|
||||||
|
|
||||||
target[0] += np.sum(np.diag(dK_dvar)*partial)
|
target[0] += np.sum(np.diag(dK_dvar)*dL_dKdiag)
|
||||||
target[1] += np.sum(np.diag(dK_dlen)*partial)
|
target[1] += np.sum(np.diag(dK_dlen)*dL_dKdiag)
|
||||||
target[2] += np.sum(np.diag(dK_dper)*partial)
|
target[2] += np.sum(np.diag(dK_dper)*dL_dKdiag)
|
||||||
|
|
|
||||||
|
|
@ -46,7 +46,7 @@ class periodic_Matern52(kernpart):
|
||||||
r = np.sqrt(r1**2 + r2**2)
|
r = np.sqrt(r1**2 + r2**2)
|
||||||
psi = np.where(r1 != 0, (np.arctan(r2/r1) + (r1<0.)*np.pi),np.arcsin(r2))
|
psi = np.where(r1 != 0, (np.arctan(r2/r1) + (r1<0.)*np.pi),np.arcsin(r2))
|
||||||
return r,omega[:,0:1], psi
|
return r,omega[:,0:1], psi
|
||||||
|
|
||||||
def _int_computation(self,r1,omega1,phi1,r2,omega2,phi2):
|
def _int_computation(self,r1,omega1,phi1,r2,omega2,phi2):
|
||||||
Gint1 = 1./(omega1+omega2.T)*( np.sin((omega1+omega2.T)*self.upper+phi1+phi2.T) - np.sin((omega1+omega2.T)*self.lower+phi1+phi2.T)) + 1./(omega1-omega2.T)*( np.sin((omega1-omega2.T)*self.upper+phi1-phi2.T) - np.sin((omega1-omega2.T)*self.lower+phi1-phi2.T) )
|
Gint1 = 1./(omega1+omega2.T)*( np.sin((omega1+omega2.T)*self.upper+phi1+phi2.T) - np.sin((omega1+omega2.T)*self.lower+phi1+phi2.T)) + 1./(omega1-omega2.T)*( np.sin((omega1-omega2.T)*self.upper+phi1-phi2.T) - np.sin((omega1-omega2.T)*self.lower+phi1-phi2.T) )
|
||||||
Gint2 = 1./(omega1+omega2.T)*( np.sin((omega1+omega2.T)*self.upper+phi1+phi2.T) - np.sin((omega1+omega2.T)*self.lower+phi1+phi2.T)) + np.cos(phi1-phi2.T)*(self.upper-self.lower)
|
Gint2 = 1./(omega1+omega2.T)*( np.sin((omega1+omega2.T)*self.upper+phi1+phi2.T) - np.sin((omega1+omega2.T)*self.lower+phi1+phi2.T)) + np.cos(phi1-phi2.T)*(self.upper-self.lower)
|
||||||
|
|
@ -105,7 +105,7 @@ class periodic_Matern52(kernpart):
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
np.add(target,np.diag(mdot(FX,self.Gi,FX.T)),target)
|
np.add(target,np.diag(mdot(FX,self.Gi,FX.T)),target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters (shape is NxMxNparam)"""
|
"""derivative of the covariance matrix with respect to the parameters (shape is NxMxNparam)"""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
|
|
@ -178,13 +178,13 @@ class periodic_Matern52(kernpart):
|
||||||
dK_dper = mdot(dFX_dper,self.Gi,FX2.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX2.T) + mdot(FX,self.Gi,dFX2_dper.T)
|
dK_dper = mdot(dFX_dper,self.Gi,FX2.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX2.T) + mdot(FX,self.Gi,dFX2_dper.T)
|
||||||
|
|
||||||
# np.add(target[:,:,0],dK_dvar, target[:,:,0])
|
# np.add(target[:,:,0],dK_dvar, target[:,:,0])
|
||||||
target[0] += np.sum(dK_dvar*partial)
|
target[0] += np.sum(dK_dvar*dL_dK)
|
||||||
#np.add(target[:,:,1],dK_dlen, target[:,:,1])
|
#np.add(target[:,:,1],dK_dlen, target[:,:,1])
|
||||||
target[1] += np.sum(dK_dlen*partial)
|
target[1] += np.sum(dK_dlen*dL_dK)
|
||||||
#np.add(target[:,:,2],dK_dper, target[:,:,2])
|
#np.add(target[:,:,2],dK_dper, target[:,:,2])
|
||||||
target[2] += np.sum(dK_dper*partial)
|
target[2] += np.sum(dK_dper*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""derivative of the diagonal of the covariance matrix with respect to the parameters"""
|
"""derivative of the diagonal of the covariance matrix with respect to the parameters"""
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
|
|
||||||
|
|
@ -251,6 +251,6 @@ class periodic_Matern52(kernpart):
|
||||||
dG_dper = 1./self.variance*(3*self.lengthscale**5/(400*np.sqrt(5))*dGint_dper + 0.5*dlower_terms_dper)
|
dG_dper = 1./self.variance*(3*self.lengthscale**5/(400*np.sqrt(5))*dGint_dper + 0.5*dlower_terms_dper)
|
||||||
dK_dper = 2*mdot(dFX_dper,self.Gi,FX.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX.T)
|
dK_dper = 2*mdot(dFX_dper,self.Gi,FX.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX.T)
|
||||||
|
|
||||||
target[0] += np.sum(np.diag(dK_dvar)*partial)
|
target[0] += np.sum(np.diag(dK_dvar)*dL_dKdiag)
|
||||||
target[1] += np.sum(np.diag(dK_dlen)*partial)
|
target[1] += np.sum(np.diag(dK_dlen)*dL_dKdiag)
|
||||||
target[2] += np.sum(np.diag(dK_dper)*partial)
|
target[2] += np.sum(np.diag(dK_dper)*dL_dKdiag)
|
||||||
|
|
|
||||||
|
|
@ -101,7 +101,7 @@ class periodic_exponential(kernpart):
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
np.add(target,np.diag(mdot(FX,self.Gi,FX.T)),target)
|
np.add(target,np.diag(mdot(FX,self.Gi,FX.T)),target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters (shape is NxMxNparam)"""
|
"""derivative of the covariance matrix with respect to the parameters (shape is NxMxNparam)"""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
|
|
@ -162,11 +162,11 @@ class periodic_exponential(kernpart):
|
||||||
|
|
||||||
dK_dper = mdot(dFX_dper,self.Gi,FX2.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX2.T) + mdot(FX,self.Gi,dFX2_dper.T)
|
dK_dper = mdot(dFX_dper,self.Gi,FX2.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX2.T) + mdot(FX,self.Gi,dFX2_dper.T)
|
||||||
|
|
||||||
target[0] += np.sum(dK_dvar*partial)
|
target[0] += np.sum(dK_dvar*dL_dK)
|
||||||
target[1] += np.sum(dK_dlen*partial)
|
target[1] += np.sum(dK_dlen*dL_dK)
|
||||||
target[2] += np.sum(dK_dper*partial)
|
target[2] += np.sum(dK_dper*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""derivative of the diagonal of the covariance matrix with respect to the parameters"""
|
"""derivative of the diagonal of the covariance matrix with respect to the parameters"""
|
||||||
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
FX = self._cos(self.basis_alpha[None,:],self.basis_omega[None,:],self.basis_phi[None,:])(X)
|
||||||
|
|
||||||
|
|
@ -222,7 +222,7 @@ class periodic_exponential(kernpart):
|
||||||
|
|
||||||
dK_dper = 2*mdot(dFX_dper,self.Gi,FX.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX.T)
|
dK_dper = 2*mdot(dFX_dper,self.Gi,FX.T) - mdot(FX,self.Gi,dG_dper,self.Gi,FX.T)
|
||||||
|
|
||||||
target[0] += np.sum(np.diag(dK_dvar)*partial)
|
target[0] += np.sum(np.diag(dK_dvar)*dL_dKdiag)
|
||||||
target[1] += np.sum(np.diag(dK_dlen)*partial)
|
target[1] += np.sum(np.diag(dK_dlen)*dL_dKdiag)
|
||||||
target[2] += np.sum(np.diag(dK_dper)*partial)
|
target[2] += np.sum(np.diag(dK_dper)*dL_dKdiag)
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -55,7 +55,7 @@ class prod(kernpart):
|
||||||
self.k2.Kdiag(X,target2)
|
self.k2.Kdiag(X,target2)
|
||||||
target += target1 * target2
|
target += target1 * target2
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters."""
|
"""derivative of the covariance matrix with respect to the parameters."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
||||||
|
|
@ -65,13 +65,13 @@ class prod(kernpart):
|
||||||
|
|
||||||
k1_target = np.zeros(self.k1.Nparam)
|
k1_target = np.zeros(self.k1.Nparam)
|
||||||
k2_target = np.zeros(self.k2.Nparam)
|
k2_target = np.zeros(self.k2.Nparam)
|
||||||
self.k1.dK_dtheta(partial*K2, X, X2, k1_target)
|
self.k1.dK_dtheta(dL_dK*K2, X, X2, k1_target)
|
||||||
self.k2.dK_dtheta(partial*K1, X, X2, k2_target)
|
self.k2.dK_dtheta(dL_dK*K1, X, X2, k2_target)
|
||||||
|
|
||||||
target[:self.k1.Nparam] += k1_target
|
target[:self.k1.Nparam] += k1_target
|
||||||
target[self.k1.Nparam:] += k2_target
|
target[self.k1.Nparam:] += k2_target
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to X."""
|
"""derivative of the covariance matrix with respect to X."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
||||||
|
|
@ -79,19 +79,19 @@ class prod(kernpart):
|
||||||
self.k1.K(X,X2,K1)
|
self.k1.K(X,X2,K1)
|
||||||
self.k2.K(X,X2,K2)
|
self.k2.K(X,X2,K2)
|
||||||
|
|
||||||
self.k1.dK_dX(partial*K2, X, X2, target)
|
self.k1.dK_dX(dL_dK*K2, X, X2, target)
|
||||||
self.k2.dK_dX(partial*K1, X, X2, target)
|
self.k2.dK_dX(dL_dK*K1, X, X2, target)
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
target1 = np.zeros((X.shape[0],))
|
target1 = np.zeros((X.shape[0],))
|
||||||
target2 = np.zeros((X.shape[0],))
|
target2 = np.zeros((X.shape[0],))
|
||||||
self.k1.Kdiag(X,target1)
|
self.k1.Kdiag(X,target1)
|
||||||
self.k2.Kdiag(X,target2)
|
self.k2.Kdiag(X,target2)
|
||||||
|
|
||||||
self.k1.dKdiag_dX(partial*target2, X, target)
|
self.k1.dKdiag_dX(dL_dKdiag*target2, X, target)
|
||||||
self.k2.dKdiag_dX(partial*target1, X, target)
|
self.k2.dKdiag_dX(dL_dKdiag*target1, X, target)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
target1 = np.zeros((X.shape[0],))
|
target1 = np.zeros((X.shape[0],))
|
||||||
target2 = np.zeros((X.shape[0],))
|
target2 = np.zeros((X.shape[0],))
|
||||||
|
|
@ -100,8 +100,8 @@ class prod(kernpart):
|
||||||
|
|
||||||
k1_target = np.zeros(self.k1.Nparam)
|
k1_target = np.zeros(self.k1.Nparam)
|
||||||
k2_target = np.zeros(self.k2.Nparam)
|
k2_target = np.zeros(self.k2.Nparam)
|
||||||
self.k1.dKdiag_dtheta(partial*target2, X, k1_target)
|
self.k1.dKdiag_dtheta(dL_dKdiag*target2, X, k1_target)
|
||||||
self.k2.dKdiag_dtheta(partial*target1, X, k2_target)
|
self.k2.dKdiag_dtheta(dL_dKdiag*target1, X, k2_target)
|
||||||
|
|
||||||
target[:self.k1.Nparam] += k1_target
|
target[:self.k1.Nparam] += k1_target
|
||||||
target[self.k1.Nparam:] += k2_target
|
target[self.k1.Nparam:] += k2_target
|
||||||
|
|
|
||||||
|
|
@ -46,7 +46,7 @@ class prod_orthogonal(kernpart):
|
||||||
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],target2)
|
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],target2)
|
||||||
target += target1 * target2
|
target += target1 * target2
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters."""
|
"""derivative of the covariance matrix with respect to the parameters."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
||||||
|
|
@ -54,8 +54,8 @@ class prod_orthogonal(kernpart):
|
||||||
self.k1.K(X[:,:self.k1.D],X2[:,:self.k1.D],K1)
|
self.k1.K(X[:,:self.k1.D],X2[:,:self.k1.D],K1)
|
||||||
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],K2)
|
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],K2)
|
||||||
|
|
||||||
self.k1.dK_dtheta(partial*K2, X[:,:self.k1.D], X2[:,:self.k1.D], target[:self.k1.Nparam])
|
self.k1.dK_dtheta(dL_dK*K2, X[:,:self.k1.D], X2[:,:self.k1.D], target[:self.k1.Nparam])
|
||||||
self.k2.dK_dtheta(partial*K1, X[:,self.k1.D:], X2[:,self.k1.D:], target[self.k1.Nparam:])
|
self.k2.dK_dtheta(dL_dK*K1, X[:,self.k1.D:], X2[:,self.k1.D:], target[self.k1.Nparam:])
|
||||||
|
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
|
|
@ -65,15 +65,15 @@ class prod_orthogonal(kernpart):
|
||||||
self.k2.Kdiag(X[:,self.k1.D:],target2)
|
self.k2.Kdiag(X[:,self.k1.D:],target2)
|
||||||
target += target1 * target2
|
target += target1 * target2
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
K1 = np.zeros(X.shape[0])
|
K1 = np.zeros(X.shape[0])
|
||||||
K2 = np.zeros(X.shape[0])
|
K2 = np.zeros(X.shape[0])
|
||||||
self.k1.Kdiag(X[:,:self.k1.D],K1)
|
self.k1.Kdiag(X[:,:self.k1.D],K1)
|
||||||
self.k2.Kdiag(X[:,self.k1.D:],K2)
|
self.k2.Kdiag(X[:,self.k1.D:],K2)
|
||||||
self.k1.dKdiag_dtheta(partial*K2,X[:,:self.k1.D],target[:self.k1.Nparam])
|
self.k1.dKdiag_dtheta(dL_dKdiag*K2,X[:,:self.k1.D],target[:self.k1.Nparam])
|
||||||
self.k2.dKdiag_dtheta(partial*K1,X[:,self.k1.D:],target[self.k1.Nparam:])
|
self.k2.dKdiag_dtheta(dL_dKdiag*K1,X[:,self.k1.D:],target[self.k1.Nparam:])
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to X."""
|
"""derivative of the covariance matrix with respect to X."""
|
||||||
if X2 is None: X2 = X
|
if X2 is None: X2 = X
|
||||||
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
K1 = np.zeros((X.shape[0],X2.shape[0]))
|
||||||
|
|
@ -81,15 +81,15 @@ class prod_orthogonal(kernpart):
|
||||||
self.k1.K(X[:,0:self.k1.D],X2[:,0:self.k1.D],K1)
|
self.k1.K(X[:,0:self.k1.D],X2[:,0:self.k1.D],K1)
|
||||||
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],K2)
|
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],K2)
|
||||||
|
|
||||||
self.k1.dK_dX(partial*K2, X[:,:self.k1.D], X2[:,:self.k1.D], target)
|
self.k1.dK_dX(dL_dK*K2, X[:,:self.k1.D], X2[:,:self.k1.D], target)
|
||||||
self.k2.dK_dX(partial*K1, X[:,self.k1.D:], X2[:,self.k1.D:], target)
|
self.k2.dK_dX(dL_dK*K1, X[:,self.k1.D:], X2[:,self.k1.D:], target)
|
||||||
|
|
||||||
def dKdiag_dX(self, partial, X, target):
|
def dKdiag_dX(self, dL_dKdiag, X, target):
|
||||||
K1 = np.zeros(X.shape[0])
|
K1 = np.zeros(X.shape[0])
|
||||||
K2 = np.zeros(X.shape[0])
|
K2 = np.zeros(X.shape[0])
|
||||||
self.k1.Kdiag(X[:,0:self.k1.D],K1)
|
self.k1.Kdiag(X[:,0:self.k1.D],K1)
|
||||||
self.k2.Kdiag(X[:,self.k1.D:],K2)
|
self.k2.Kdiag(X[:,self.k1.D:],K2)
|
||||||
|
|
||||||
self.k1.dK_dX(partial*K2, X[:,:self.k1.D], target)
|
self.k1.dK_dX(dL_dKdiag*K2, X[:,:self.k1.D], target)
|
||||||
self.k2.dK_dX(partial*K1, X[:,self.k1.D:], target)
|
self.k2.dK_dX(dL_dKdiag*K1, X[:,self.k1.D:], target)
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -82,27 +82,27 @@ class rbf(kernpart):
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
np.add(target,self.variance,target)
|
np.add(target,self.variance,target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
self._K_computations(X,X2)
|
self._K_computations(X,X2)
|
||||||
target[0] += np.sum(self._K_dvar*partial)
|
target[0] += np.sum(self._K_dvar*dL_dK)
|
||||||
if self.ARD == True:
|
if self.ARD == True:
|
||||||
dl = self._K_dvar[:,:,None]*self.variance*self._K_dist2/self.lengthscale
|
dl = self._K_dvar[:,:,None]*self.variance*self._K_dist2/self.lengthscale
|
||||||
target[1:] += (dl*partial[:,:,None]).sum(0).sum(0)
|
target[1:] += (dl*dL_dK[:,:,None]).sum(0).sum(0)
|
||||||
else:
|
else:
|
||||||
target[1] += np.sum(self._K_dvar*self.variance*(self._K_dist2.sum(-1))/self.lengthscale*partial)
|
target[1] += np.sum(self._K_dvar*self.variance*(self._K_dist2.sum(-1))/self.lengthscale*dL_dK)
|
||||||
#np.sum(self._K_dvar*self.variance*self._K_dist2/self.lengthscale*partial)
|
#np.sum(self._K_dvar*self.variance*self._K_dist2/self.lengthscale*dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
#NB: derivative of diagonal elements wrt lengthscale is 0
|
#NB: derivative of diagonal elements wrt lengthscale is 0
|
||||||
target[0] += np.sum(partial)
|
target[0] += np.sum(dL_dKdiag)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
self._K_computations(X,X2)
|
self._K_computations(X,X2)
|
||||||
_K_dist = X[:,None,:]-X2[None,:,:]
|
_K_dist = X[:,None,:]-X2[None,:,:]
|
||||||
dK_dX = np.transpose(-self.variance*self._K_dvar[:,:,np.newaxis]*_K_dist/self.lengthscale2,(1,0,2))
|
dK_dX = np.transpose(-self.variance*self._K_dvar[:,:,np.newaxis]*_K_dist/self.lengthscale2,(1,0,2))
|
||||||
target += np.sum(dK_dX*partial.T[:,:,None],0)
|
target += np.sum(dK_dX*dL_dK.T[:,:,None],0)
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
|
|
||||||
|
|
@ -113,69 +113,69 @@ class rbf(kernpart):
|
||||||
def psi0(self,Z,mu,S,target):
|
def psi0(self,Z,mu,S,target):
|
||||||
target += self.variance
|
target += self.variance
|
||||||
|
|
||||||
def dpsi0_dtheta(self,partial,Z,mu,S,target):
|
def dpsi0_dtheta(self,dL_dpsi0,Z,mu,S,target):
|
||||||
target[0] += np.sum(partial)
|
target[0] += np.sum(dL_dpsi0)
|
||||||
|
|
||||||
def dpsi0_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi0_dmuS(self,dL_dpsi0,Z,mu,S,target_mu,target_S):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def psi1(self,Z,mu,S,target):
|
def psi1(self,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
target += self._psi1
|
target += self._psi1
|
||||||
|
|
||||||
def dpsi1_dtheta(self,partial,Z,mu,S,target):
|
def dpsi1_dtheta(self,dL_dpsi1,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
denom_deriv = S[:,None,:]/(self.lengthscale**3+self.lengthscale*S[:,None,:])
|
denom_deriv = S[:,None,:]/(self.lengthscale**3+self.lengthscale*S[:,None,:])
|
||||||
d_length = self._psi1[:,:,None]*(self.lengthscale*np.square(self._psi1_dist/(self.lengthscale2+S[:,None,:])) + denom_deriv)
|
d_length = self._psi1[:,:,None]*(self.lengthscale*np.square(self._psi1_dist/(self.lengthscale2+S[:,None,:])) + denom_deriv)
|
||||||
target[0] += np.sum(partial*self._psi1/self.variance)
|
target[0] += np.sum(dL_dpsi1*self._psi1/self.variance)
|
||||||
dpsi1_dlength = d_length*partial[:,:,None]
|
dpsi1_dlength = d_length*dL_dpsi1[:,:,None]
|
||||||
if not self.ARD:
|
if not self.ARD:
|
||||||
target[1] += dpsi1_dlength.sum()
|
target[1] += dpsi1_dlength.sum()
|
||||||
else:
|
else:
|
||||||
target[1:] += dpsi1_dlength.sum(0).sum(0)
|
target[1:] += dpsi1_dlength.sum(0).sum(0)
|
||||||
|
|
||||||
def dpsi1_dZ(self,partial,Z,mu,S,target):
|
def dpsi1_dZ(self,dL_dpsi1,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
denominator = (self.lengthscale2*(self._psi1_denom))
|
denominator = (self.lengthscale2*(self._psi1_denom))
|
||||||
dpsi1_dZ = - self._psi1[:,:,None] * ((self._psi1_dist/denominator))
|
dpsi1_dZ = - self._psi1[:,:,None] * ((self._psi1_dist/denominator))
|
||||||
target += np.sum(partial.T[:,:,None] * dpsi1_dZ, 0)
|
target += np.sum(dL_dpsi1.T[:,:,None] * dpsi1_dZ, 0)
|
||||||
|
|
||||||
def dpsi1_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi1_dmuS(self,dL_dpsi1,Z,mu,S,target_mu,target_S):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
tmp = self._psi1[:,:,None]/self.lengthscale2/self._psi1_denom
|
tmp = self._psi1[:,:,None]/self.lengthscale2/self._psi1_denom
|
||||||
target_mu += np.sum(partial.T[:, :, None]*tmp*self._psi1_dist,1)
|
target_mu += np.sum(dL_dpsi1.T[:, :, None]*tmp*self._psi1_dist,1)
|
||||||
target_S += np.sum(partial.T[:, :, None]*0.5*tmp*(self._psi1_dist_sq-1),1)
|
target_S += np.sum(dL_dpsi1.T[:, :, None]*0.5*tmp*(self._psi1_dist_sq-1),1)
|
||||||
|
|
||||||
def psi2(self,Z,mu,S,target):
|
def psi2(self,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
target += self._psi2
|
target += self._psi2
|
||||||
|
|
||||||
def dpsi2_dtheta(self,partial,Z,mu,S,target):
|
def dpsi2_dtheta(self,dL_dpsi2,Z,mu,S,target):
|
||||||
"""Shape N,M,M,Ntheta"""
|
"""Shape N,M,M,Ntheta"""
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
d_var = 2.*self._psi2/self.variance
|
d_var = 2.*self._psi2/self.variance
|
||||||
d_length = self._psi2[:,:,:,None]*(0.5*self._psi2_Zdist_sq*self._psi2_denom + 2.*self._psi2_mudist_sq + 2.*S[:,None,None,:]/self.lengthscale2)/(self.lengthscale*self._psi2_denom)
|
d_length = self._psi2[:,:,:,None]*(0.5*self._psi2_Zdist_sq*self._psi2_denom + 2.*self._psi2_mudist_sq + 2.*S[:,None,None,:]/self.lengthscale2)/(self.lengthscale*self._psi2_denom)
|
||||||
|
|
||||||
target[0] += np.sum(partial*d_var)
|
target[0] += np.sum(dL_dpsi2*d_var)
|
||||||
dpsi2_dlength = d_length*partial[:,:,:,None]
|
dpsi2_dlength = d_length*dL_dpsi2[:,:,:,None]
|
||||||
if not self.ARD:
|
if not self.ARD:
|
||||||
target[1] += dpsi2_dlength.sum()
|
target[1] += dpsi2_dlength.sum()
|
||||||
else:
|
else:
|
||||||
target[1:] += dpsi2_dlength.sum(0).sum(0).sum(0)
|
target[1:] += dpsi2_dlength.sum(0).sum(0).sum(0)
|
||||||
|
|
||||||
def dpsi2_dZ(self,partial,Z,mu,S,target):
|
def dpsi2_dZ(self,dL_dpsi2,Z,mu,S,target):
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
term1 = 0.5*self._psi2_Zdist/self.lengthscale2 # M, M, Q
|
term1 = 0.5*self._psi2_Zdist/self.lengthscale2 # M, M, Q
|
||||||
term2 = self._psi2_mudist/self._psi2_denom/self.lengthscale2 # N, M, M, Q
|
term2 = self._psi2_mudist/self._psi2_denom/self.lengthscale2 # N, M, M, Q
|
||||||
dZ = self._psi2[:,:,:,None] * (term1[None] + term2)
|
dZ = self._psi2[:,:,:,None] * (term1[None] + term2)
|
||||||
target += (partial[:,:,:,None]*dZ).sum(0).sum(0)
|
target += (dL_dpsi2[:,:,:,None]*dZ).sum(0).sum(0)
|
||||||
|
|
||||||
def dpsi2_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi2_dmuS(self,dL_dpsi2,Z,mu,S,target_mu,target_S):
|
||||||
"""Think N,M,M,Q """
|
"""Think N,M,M,Q """
|
||||||
self._psi_computations(Z,mu,S)
|
self._psi_computations(Z,mu,S)
|
||||||
tmp = self._psi2[:,:,:,None]/self.lengthscale2/self._psi2_denom
|
tmp = self._psi2[:,:,:,None]/self.lengthscale2/self._psi2_denom
|
||||||
target_mu += (partial[:,:,:,None]*-tmp*2.*self._psi2_mudist).sum(1).sum(1)
|
target_mu += (dL_dpsi2[:,:,:,None]*-tmp*2.*self._psi2_mudist).sum(1).sum(1)
|
||||||
target_S += (partial[:,:,:,None]*tmp*(2.*self._psi2_mudist_sq-1)).sum(1).sum(1)
|
target_S += (dL_dpsi2[:,:,:,None]*tmp*(2.*self._psi2_mudist_sq-1)).sum(1).sum(1)
|
||||||
|
|
||||||
|
|
||||||
#---------------------------------------#
|
#---------------------------------------#
|
||||||
|
|
|
||||||
|
|
@ -51,7 +51,7 @@ class symmetric(kernpart):
|
||||||
self.k.K(X,AX2,target)
|
self.k.K(X,AX2,target)
|
||||||
self.k.K(AX,AX2,target)
|
self.k.K(AX,AX2,target)
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to the parameters."""
|
"""derivative of the covariance matrix with respect to the parameters."""
|
||||||
AX = np.dot(X,self.transform)
|
AX = np.dot(X,self.transform)
|
||||||
if X2 is None:
|
if X2 is None:
|
||||||
|
|
@ -59,13 +59,13 @@ class symmetric(kernpart):
|
||||||
ZX2 = AX
|
ZX2 = AX
|
||||||
else:
|
else:
|
||||||
AX2 = np.dot(X2, self.transform)
|
AX2 = np.dot(X2, self.transform)
|
||||||
self.k.dK_dtheta(partial,X,X2,target)
|
self.k.dK_dtheta(dL_dK,X,X2,target)
|
||||||
self.k.dK_dtheta(partial,AX,X2,target)
|
self.k.dK_dtheta(dL_dK,AX,X2,target)
|
||||||
self.k.dK_dtheta(partial,X,AX2,target)
|
self.k.dK_dtheta(dL_dK,X,AX2,target)
|
||||||
self.k.dK_dtheta(partial,AX,AX2,target)
|
self.k.dK_dtheta(dL_dK,AX,AX2,target)
|
||||||
|
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
"""derivative of the covariance matrix with respect to X."""
|
"""derivative of the covariance matrix with respect to X."""
|
||||||
AX = np.dot(X,self.transform)
|
AX = np.dot(X,self.transform)
|
||||||
if X2 is None:
|
if X2 is None:
|
||||||
|
|
@ -73,10 +73,10 @@ class symmetric(kernpart):
|
||||||
ZX2 = AX
|
ZX2 = AX
|
||||||
else:
|
else:
|
||||||
AX2 = np.dot(X2, self.transform)
|
AX2 = np.dot(X2, self.transform)
|
||||||
self.k.dK_dX(partial, X, X2, target)
|
self.k.dK_dX(dL_dK, X, X2, target)
|
||||||
self.k.dK_dX(partial, AX, X2, target)
|
self.k.dK_dX(dL_dK, AX, X2, target)
|
||||||
self.k.dK_dX(partial, X, AX2, target)
|
self.k.dK_dX(dL_dK, X, AX2, target)
|
||||||
self.k.dK_dX(partial, AX ,AX2, target)
|
self.k.dK_dX(dL_dK, AX ,AX2, target)
|
||||||
|
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
|
|
@ -84,9 +84,9 @@ class symmetric(kernpart):
|
||||||
self.K(X,X,foo)
|
self.K(X,X,foo)
|
||||||
target += np.diag(foo)
|
target += np.diag(foo)
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
"""Compute the diagonal of the covariance matrix associated to X."""
|
"""Compute the diagonal of the covariance matrix associated to X."""
|
||||||
raise NotImplementedError
|
raise NotImplementedError
|
||||||
|
|
|
||||||
|
|
@ -37,50 +37,50 @@ class white(kernpart):
|
||||||
def Kdiag(self,X,target):
|
def Kdiag(self,X,target):
|
||||||
target += self.variance
|
target += self.variance
|
||||||
|
|
||||||
def dK_dtheta(self,partial,X,X2,target):
|
def dK_dtheta(self,dL_dK,X,X2,target):
|
||||||
if X.shape==X2.shape:
|
if X.shape==X2.shape:
|
||||||
if np.all(X==X2):
|
if np.all(X==X2):
|
||||||
target += np.trace(partial)
|
target += np.trace(dL_dK)
|
||||||
|
|
||||||
def dKdiag_dtheta(self,partial,X,target):
|
def dKdiag_dtheta(self,dL_dKdiag,X,target):
|
||||||
target += np.sum(partial)
|
target += np.sum(dL_dKdiag)
|
||||||
|
|
||||||
def dK_dX(self,partial,X,X2,target):
|
def dK_dX(self,dL_dK,X,X2,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dKdiag_dX(self,partial,X,target):
|
def dKdiag_dX(self,dL_dKdiag,X,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def psi0(self,Z,mu,S,target):
|
def psi0(self,Z,mu,S,target):
|
||||||
target += self.variance
|
target += self.variance
|
||||||
|
|
||||||
def dpsi0_dtheta(self,partial,Z,mu,S,target):
|
def dpsi0_dtheta(self,dL_dpsi0,Z,mu,S,target):
|
||||||
target += partial.sum()
|
target += dL_dpsi0.sum()
|
||||||
|
|
||||||
def dpsi0_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi0_dmuS(self,dL_dpsi0,Z,mu,S,target_mu,target_S):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def psi1(self,Z,mu,S,target):
|
def psi1(self,Z,mu,S,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi1_dtheta(self,partial,Z,mu,S,target):
|
def dpsi1_dtheta(self,dL_dpsi1,Z,mu,S,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi1_dZ(self,partial,Z,mu,S,target):
|
def dpsi1_dZ(self,dL_dpsi1,Z,mu,S,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi1_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi1_dmuS(self,dL_dpsi1,Z,mu,S,target_mu,target_S):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def psi2(self,Z,mu,S,target):
|
def psi2(self,Z,mu,S,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi2_dZ(self,partial,Z,mu,S,target):
|
def dpsi2_dZ(self,dL_dpsi2,Z,mu,S,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi2_dtheta(self,partial,Z,mu,S,target):
|
def dpsi2_dtheta(self,dL_dpsi2,Z,mu,S,target):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
def dpsi2_dmuS(self,partial,Z,mu,S,target_mu,target_S):
|
def dpsi2_dmuS(self,dL_dpsi2,Z,mu,S,target_mu,target_S):
|
||||||
pass
|
pass
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -110,7 +110,7 @@ class EP(likelihood):
|
||||||
Sroot_tilde_K = np.sqrt(self.tau_tilde)[:,None]*K
|
Sroot_tilde_K = np.sqrt(self.tau_tilde)[:,None]*K
|
||||||
B = np.eye(self.N) + np.sqrt(self.tau_tilde)[None,:]*Sroot_tilde_K
|
B = np.eye(self.N) + np.sqrt(self.tau_tilde)[None,:]*Sroot_tilde_K
|
||||||
L = jitchol(B)
|
L = jitchol(B)
|
||||||
V,info = linalg.flapack.dtrtrs(L,Sroot_tilde_K,lower=1)
|
V,info = linalg.lapack.flapack.dtrtrs(L,Sroot_tilde_K,lower=1)
|
||||||
Sigma = K - np.dot(V.T,V)
|
Sigma = K - np.dot(V.T,V)
|
||||||
mu = np.dot(Sigma,self.v_tilde)
|
mu = np.dot(Sigma,self.v_tilde)
|
||||||
epsilon_np1 = sum((self.tau_tilde-self.np1[-1])**2)/self.N
|
epsilon_np1 = sum((self.tau_tilde-self.np1[-1])**2)/self.N
|
||||||
|
|
@ -187,7 +187,7 @@ class EP(likelihood):
|
||||||
#Posterior distribution parameters update
|
#Posterior distribution parameters update
|
||||||
LLT = LLT + np.outer(Kmn[:,i],Kmn[:,i])*Delta_tau
|
LLT = LLT + np.outer(Kmn[:,i],Kmn[:,i])*Delta_tau
|
||||||
L = jitchol(LLT)
|
L = jitchol(LLT)
|
||||||
V,info = linalg.flapack.dtrtrs(L,Kmn,lower=1)
|
V,info = linalg.lapack.flapack.dtrtrs(L,Kmn,lower=1)
|
||||||
Sigma_diag = np.sum(V*V,-2)
|
Sigma_diag = np.sum(V*V,-2)
|
||||||
si = np.sum(V.T*V[:,i],-1)
|
si = np.sum(V.T*V[:,i],-1)
|
||||||
mu = mu + (Delta_v-Delta_tau*mu[i])*si
|
mu = mu + (Delta_v-Delta_tau*mu[i])*si
|
||||||
|
|
@ -195,8 +195,8 @@ class EP(likelihood):
|
||||||
#Sigma recomputation with Cholesky decompositon
|
#Sigma recomputation with Cholesky decompositon
|
||||||
LLT0 = LLT0 + np.dot(Kmn*self.tau_tilde[None,:],Kmn.T)
|
LLT0 = LLT0 + np.dot(Kmn*self.tau_tilde[None,:],Kmn.T)
|
||||||
L = jitchol(LLT)
|
L = jitchol(LLT)
|
||||||
V,info = linalg.flapack.dtrtrs(L,Kmn,lower=1)
|
V,info = linalg.lapack.flapack.dtrtrs(L,Kmn,lower=1)
|
||||||
V2,info = linalg.flapack.dtrtrs(L.T,V,lower=0)
|
V2,info = linalg.lapack.flapack.dtrtrs(L.T,V,lower=0)
|
||||||
Sigma_diag = np.sum(V*V,-2)
|
Sigma_diag = np.sum(V*V,-2)
|
||||||
Knmv_tilde = np.dot(Kmn,self.v_tilde)
|
Knmv_tilde = np.dot(Kmn,self.v_tilde)
|
||||||
mu = np.dot(V2.T,Knmv_tilde)
|
mu = np.dot(V2.T,Knmv_tilde)
|
||||||
|
|
@ -295,7 +295,7 @@ class EP(likelihood):
|
||||||
P = (Diag / Diag0)[:,None] * P0
|
P = (Diag / Diag0)[:,None] * P0
|
||||||
RPT0 = np.dot(R0,P0.T)
|
RPT0 = np.dot(R0,P0.T)
|
||||||
L = jitchol(np.eye(M) + np.dot(RPT0,(1./Diag0 - Diag/(Diag0**2))[:,None]*RPT0.T))
|
L = jitchol(np.eye(M) + np.dot(RPT0,(1./Diag0 - Diag/(Diag0**2))[:,None]*RPT0.T))
|
||||||
R,info = linalg.flapack.dtrtrs(L,R0,lower=1)
|
R,info = linalg.lapack.flapack.dtrtrs(L,R0,lower=1)
|
||||||
RPT = np.dot(R,P.T)
|
RPT = np.dot(R,P.T)
|
||||||
Sigma_diag = Diag + np.sum(RPT.T*RPT.T,-1)
|
Sigma_diag = Diag + np.sum(RPT.T*RPT.T,-1)
|
||||||
self.w = Diag * self.v_tilde
|
self.w = Diag * self.v_tilde
|
||||||
|
|
|
||||||
|
|
@ -8,7 +8,7 @@ class Gaussian(likelihood):
|
||||||
self.Z = 0. # a correction factor which accounts for the approximation made
|
self.Z = 0. # a correction factor which accounts for the approximation made
|
||||||
N, self.D = data.shape
|
N, self.D = data.shape
|
||||||
|
|
||||||
#normalisation
|
#normaliztion
|
||||||
if normalize:
|
if normalize:
|
||||||
self._mean = data.mean(0)[None,:]
|
self._mean = data.mean(0)[None,:]
|
||||||
self._std = data.std(0)[None,:]
|
self._std = data.std(0)[None,:]
|
||||||
|
|
@ -45,7 +45,7 @@ class Gaussian(likelihood):
|
||||||
|
|
||||||
def predictive_values(self,mu,var):
|
def predictive_values(self,mu,var):
|
||||||
"""
|
"""
|
||||||
Un-normalise the prediction and add the likelihood variance, then return the 5%, 95% interval
|
Un-normalize the prediction and add the likelihood variance, then return the 5%, 95% interval
|
||||||
"""
|
"""
|
||||||
mean = mu*self._std + self._mean
|
mean = mu*self._std + self._mean
|
||||||
true_var = (var + self._variance)*self._std**2
|
true_var = (var + self._variance)*self._std**2
|
||||||
|
|
|
||||||
|
|
@ -30,7 +30,6 @@ class GP(model):
|
||||||
.. Note:: Multiple independent outputs are allowed using columns of Y
|
.. Note:: Multiple independent outputs are allowed using columns of Y
|
||||||
|
|
||||||
"""
|
"""
|
||||||
#FIXME normalize vs normalise
|
|
||||||
def __init__(self, X, likelihood, kernel, normalize_X=False, Xslices=None):
|
def __init__(self, X, likelihood, kernel, normalize_X=False, Xslices=None):
|
||||||
|
|
||||||
# parse arguments
|
# parse arguments
|
||||||
|
|
@ -41,7 +40,7 @@ class GP(model):
|
||||||
assert isinstance(kernel, kern.kern)
|
assert isinstance(kernel, kern.kern)
|
||||||
self.kern = kernel
|
self.kern = kernel
|
||||||
|
|
||||||
#here's some simple normalisation for the inputs
|
#here's some simple normalization for the inputs
|
||||||
if normalize_X:
|
if normalize_X:
|
||||||
self._Xmean = X.mean(0)[None,:]
|
self._Xmean = X.mean(0)[None,:]
|
||||||
self._Xstd = X.std(0)[None,:]
|
self._Xstd = X.std(0)[None,:]
|
||||||
|
|
@ -129,12 +128,12 @@ class GP(model):
|
||||||
|
|
||||||
For the likelihood parameters, pass in alpha = K^-1 y
|
For the likelihood parameters, pass in alpha = K^-1 y
|
||||||
"""
|
"""
|
||||||
return np.hstack((self.kern.dK_dtheta(partial=self.dL_dK,X=self.X,slices1=self.Xslices,slices2=self.Xslices), self.likelihood._gradients(partial=np.diag(self.dL_dK))))
|
return np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK,X=self.X,slices1=self.Xslices,slices2=self.Xslices), self.likelihood._gradients(partial=np.diag(self.dL_dK))))
|
||||||
|
|
||||||
def _raw_predict(self,_Xnew,slices=None, full_cov=False):
|
def _raw_predict(self,_Xnew,slices=None, full_cov=False):
|
||||||
"""
|
"""
|
||||||
Internal helper function for making predictions, does not account
|
Internal helper function for making predictions, does not account
|
||||||
for normalisation or likelihood
|
for normalization or likelihood
|
||||||
"""
|
"""
|
||||||
Kx = self.kern.K(self.X,_Xnew, slices1=self.Xslices,slices2=slices)
|
Kx = self.kern.K(self.X,_Xnew, slices1=self.Xslices,slices2=slices)
|
||||||
mu = np.dot(np.dot(Kx.T,self.Ki),self.likelihood.Y)
|
mu = np.dot(np.dot(Kx.T,self.Ki),self.likelihood.Y)
|
||||||
|
|
@ -172,10 +171,10 @@ class GP(model):
|
||||||
- If a list of booleans, specifying which kernel parts are active
|
- If a list of booleans, specifying which kernel parts are active
|
||||||
|
|
||||||
If full_cov and self.D > 1, the return shape of var is Nnew x Nnew x self.D. If self.D == 1, the return shape is Nnew x Nnew.
|
If full_cov and self.D > 1, the return shape of var is Nnew x Nnew x self.D. If self.D == 1, the return shape is Nnew x Nnew.
|
||||||
This is to allow for different normalisations of the output dimensions.
|
This is to allow for different normalizations of the output dimensions.
|
||||||
|
|
||||||
"""
|
"""
|
||||||
#normalise X values
|
#normalize X values
|
||||||
Xnew = (Xnew.copy() - self._Xmean) / self._Xstd
|
Xnew = (Xnew.copy() - self._Xmean) / self._Xstd
|
||||||
mu, var = self._raw_predict(Xnew, slices, full_cov)
|
mu, var = self._raw_predict(Xnew, slices, full_cov)
|
||||||
|
|
||||||
|
|
@ -187,7 +186,7 @@ class GP(model):
|
||||||
|
|
||||||
def plot_f(self, samples=0, plot_limits=None, which_data='all', which_functions='all', resolution=None, full_cov=False):
|
def plot_f(self, samples=0, plot_limits=None, which_data='all', which_functions='all', resolution=None, full_cov=False):
|
||||||
"""
|
"""
|
||||||
Plot the GP's view of the world, where the data is normalised and the likelihood is Gaussian
|
Plot the GP's view of the world, where the data is normalized and the likelihood is Gaussian
|
||||||
|
|
||||||
:param samples: the number of a posteriori samples to plot
|
:param samples: the number of a posteriori samples to plot
|
||||||
:param which_data: which if the training data to plot (default all)
|
:param which_data: which if the training data to plot (default all)
|
||||||
|
|
@ -203,7 +202,7 @@ class GP(model):
|
||||||
- In higher dimensions, we've no implemented this yet !TODO!
|
- In higher dimensions, we've no implemented this yet !TODO!
|
||||||
|
|
||||||
Can plot only part of the data and part of the posterior functions using which_data and which_functions
|
Can plot only part of the data and part of the posterior functions using which_data and which_functions
|
||||||
Plot the data's view of the world, with non-normalised values and GP predictions passed through the likelihood
|
Plot the data's view of the world, with non-normalized values and GP predictions passed through the likelihood
|
||||||
"""
|
"""
|
||||||
if which_functions=='all':
|
if which_functions=='all':
|
||||||
which_functions = [True]*self.kern.Nparts
|
which_functions = [True]*self.kern.Nparts
|
||||||
|
|
@ -221,7 +220,7 @@ class GP(model):
|
||||||
Ysim = np.random.multivariate_normal(m.flatten(),v,samples)
|
Ysim = np.random.multivariate_normal(m.flatten(),v,samples)
|
||||||
gpplot(Xnew,m,m-2*np.sqrt(np.diag(v)[:,None]),m+2*np.sqrt(np.diag(v))[:,None])
|
gpplot(Xnew,m,m-2*np.sqrt(np.diag(v)[:,None]),m+2*np.sqrt(np.diag(v))[:,None])
|
||||||
for i in range(samples):
|
for i in range(samples):
|
||||||
pb.plot(Xnew,Ysim[i,:],Tango.coloursHex['darkBlue'],linewidth=0.25)
|
pb.plot(Xnew,Ysim[i,:],Tango.colorsHex['darkBlue'],linewidth=0.25)
|
||||||
pb.plot(self.X[which_data],self.likelihood.Y[which_data],'kx',mew=1.5)
|
pb.plot(self.X[which_data],self.likelihood.Y[which_data],'kx',mew=1.5)
|
||||||
pb.xlim(xmin,xmax)
|
pb.xlim(xmin,xmax)
|
||||||
ymin,ymax = min(np.append(self.likelihood.Y,m-2*np.sqrt(np.diag(v)[:,None]))), max(np.append(self.likelihood.Y,m+2*np.sqrt(np.diag(v)[:,None])))
|
ymin,ymax = min(np.append(self.likelihood.Y,m-2*np.sqrt(np.diag(v)[:,None]))), max(np.append(self.likelihood.Y,m+2*np.sqrt(np.diag(v)[:,None])))
|
||||||
|
|
|
||||||
|
|
@ -54,7 +54,7 @@ class sparse_GP(GP):
|
||||||
|
|
||||||
GP.__init__(self, X, likelihood, kernel=kernel, normalize_X=normalize_X, Xslices=Xslices)
|
GP.__init__(self, X, likelihood, kernel=kernel, normalize_X=normalize_X, Xslices=Xslices)
|
||||||
|
|
||||||
#normalise X uncertainty also
|
#normalize X uncertainty also
|
||||||
if self.has_uncertain_inputs:
|
if self.has_uncertain_inputs:
|
||||||
self.X_uncertainty /= np.square(self._Xstd)
|
self.X_uncertainty /= np.square(self._Xstd)
|
||||||
|
|
||||||
|
|
@ -208,7 +208,7 @@ class sparse_GP(GP):
|
||||||
if self.has_uncertain_inputs:
|
if self.has_uncertain_inputs:
|
||||||
dL_dtheta += self.kern.dpsi0_dtheta(self.dL_dpsi0, self.Z,self.X,self.X_uncertainty)
|
dL_dtheta += self.kern.dpsi0_dtheta(self.dL_dpsi0, self.Z,self.X,self.X_uncertainty)
|
||||||
dL_dtheta += self.kern.dpsi1_dtheta(self.dL_dpsi1.T,self.Z,self.X, self.X_uncertainty)
|
dL_dtheta += self.kern.dpsi1_dtheta(self.dL_dpsi1.T,self.Z,self.X, self.X_uncertainty)
|
||||||
dL_dtheta += self.kern.dpsi2_dtheta(self.dL_dpsi2,self.dL_dpsi1.T, self.Z,self.X, self.X_uncertainty)
|
dL_dtheta += self.kern.dpsi2_dtheta(self.dL_dpsi2, self.Z,self.X, self.X_uncertainty)
|
||||||
else:
|
else:
|
||||||
dL_dtheta += self.kern.dK_dtheta(self.dL_dpsi1,self.Z,self.X)
|
dL_dtheta += self.kern.dK_dtheta(self.dL_dpsi1,self.Z,self.X)
|
||||||
dL_dtheta += self.kern.dKdiag_dtheta(self.dL_dpsi0, self.X)
|
dL_dtheta += self.kern.dKdiag_dtheta(self.dL_dpsi0, self.X)
|
||||||
|
|
@ -228,7 +228,7 @@ class sparse_GP(GP):
|
||||||
return dL_dZ
|
return dL_dZ
|
||||||
|
|
||||||
def _raw_predict(self, Xnew, slices, full_cov=False):
|
def _raw_predict(self, Xnew, slices, full_cov=False):
|
||||||
"""Internal helper function for making predictions, does not account for normalisation"""
|
"""Internal helper function for making predictions, does not account for normalization"""
|
||||||
|
|
||||||
Kx = self.kern.K(self.Z, Xnew)
|
Kx = self.kern.K(self.Z, Xnew)
|
||||||
mu = mdot(Kx.T, self.C/self.scale_factor, self.psi1V)
|
mu = mdot(Kx.T, self.C/self.scale_factor, self.psi1V)
|
||||||
|
|
|
||||||
|
|
@ -43,7 +43,7 @@ class sparse_GPLVM(sparse_GP_regression, GPLVM):
|
||||||
|
|
||||||
def dL_dX(self):
|
def dL_dX(self):
|
||||||
dL_dX = self.kern.dKdiag_dX(self.dL_dpsi0,self.X)
|
dL_dX = self.kern.dKdiag_dX(self.dL_dpsi0,self.X)
|
||||||
dL_dX += self.kern.dK_dX(self.dL_dpsi1,self.X,self.Z)
|
dL_dX += self.kern.dK_dX(self.dL_dpsi1.T,self.X,self.Z)
|
||||||
|
|
||||||
return dL_dX
|
return dL_dX
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -93,7 +93,7 @@ class uncollapsed_sparse_GP(sparse_GP):
|
||||||
return A+B+C+D+E
|
return A+B+C+D+E
|
||||||
|
|
||||||
def _raw_predict(self, Xnew, slices,full_cov=False):
|
def _raw_predict(self, Xnew, slices,full_cov=False):
|
||||||
"""Internal helper function for making predictions, does not account for normalisation"""
|
"""Internal helper function for making predictions, does not account for normalization"""
|
||||||
Kx = self.kern.K(Xnew,self.Z)
|
Kx = self.kern.K(Xnew,self.Z)
|
||||||
mu = mdot(Kx,self.Kmmi,self.q_u_expectation[0])
|
mu = mdot(Kx,self.Kmmi,self.q_u_expectation[0])
|
||||||
|
|
||||||
|
|
|
||||||
0
GPy/testing/__init__.py
Normal file
0
GPy/testing/__init__.py
Normal file
|
|
@ -12,6 +12,7 @@ class BGPLVMTests(unittest.TestCase):
|
||||||
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
K = k.K(X)
|
K = k.K(X)
|
||||||
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
Y -= Y.mean(axis=0)
|
||||||
k = GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
k = GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
m.constrain_positive('(rbf|bias|noise|white|S)')
|
m.constrain_positive('(rbf|bias|noise|white|S)')
|
||||||
|
|
@ -24,6 +25,7 @@ class BGPLVMTests(unittest.TestCase):
|
||||||
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
K = k.K(X)
|
K = k.K(X)
|
||||||
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
Y -= Y.mean(axis=0)
|
||||||
k = GPy.kern.linear(Q) + GPy.kern.white(Q, 0.00001)
|
k = GPy.kern.linear(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
m.constrain_positive('(linear|bias|noise|white|S)')
|
m.constrain_positive('(linear|bias|noise|white|S)')
|
||||||
|
|
@ -36,13 +38,40 @@ class BGPLVMTests(unittest.TestCase):
|
||||||
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
K = k.K(X)
|
K = k.K(X)
|
||||||
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
Y -= Y.mean(axis=0)
|
||||||
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
m.constrain_positive('(rbf|bias|noise|white|S)')
|
m.constrain_positive('(rbf|bias|noise|white|S)')
|
||||||
m.randomize()
|
m.randomize()
|
||||||
self.assertTrue(m.checkgrad())
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
def test_rbf_bias_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
Y -= Y.mean(axis=0)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
|
m.constrain_positive('(rbf|bias|noise|white|S)')
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
def test_linear_bias_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.linear(Q) + GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
Y -= Y.mean(axis=0)
|
||||||
|
k = GPy.kern.linear(Q) + GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.Bayesian_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
|
m.constrain_positive('(linear|bias|noise|white|S)')
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
|
||||||
if __name__ == "__main__":
|
if __name__ == "__main__":
|
||||||
print "Running unit tests, please be (very) patient..."
|
print "Running unit tests, please be (very) patient..."
|
||||||
unittest.main()
|
unittest.main()
|
||||||
|
|
|
||||||
26
GPy/testing/examples_tests.py
Normal file
26
GPy/testing/examples_tests.py
Normal file
|
|
@ -0,0 +1,26 @@
|
||||||
|
# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
|
||||||
|
# Licensed under the BSD 3-clause license (see LICENSE.txt)
|
||||||
|
|
||||||
|
import unittest
|
||||||
|
import numpy as np
|
||||||
|
import GPy
|
||||||
|
|
||||||
|
class ExamplesTests(unittest.TestCase):
|
||||||
|
def test_check_model_returned(self):
|
||||||
|
pass
|
||||||
|
|
||||||
|
def test_model_checkgrads(self):
|
||||||
|
pass
|
||||||
|
|
||||||
|
def test_all_examples(self):
|
||||||
|
pass
|
||||||
|
#Load models
|
||||||
|
|
||||||
|
#Loop through models
|
||||||
|
#for model in models:
|
||||||
|
#self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
|
||||||
|
if __name__ == "__main__":
|
||||||
|
print "Running unit tests, please be (very) patient..."
|
||||||
|
unittest.main()
|
||||||
47
GPy/testing/gplvm_tests.py
Normal file
47
GPy/testing/gplvm_tests.py
Normal file
|
|
@ -0,0 +1,47 @@
|
||||||
|
# Copyright (c) 2012, Nicolo Fusi
|
||||||
|
# Licensed under the BSD 3-clause license (see LICENSE.txt)
|
||||||
|
|
||||||
|
import unittest
|
||||||
|
import numpy as np
|
||||||
|
import GPy
|
||||||
|
|
||||||
|
class GPLVMTests(unittest.TestCase):
|
||||||
|
def test_bias_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
k = GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.GPLVM(Y, Q, kernel = k)
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
def test_linear_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
k = GPy.kern.linear(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.GPLVM(Y, Q, kernel = k)
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
def test_rbf_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.GPLVM(Y, Q, kernel = k)
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
if __name__ == "__main__":
|
||||||
|
print "Running unit tests, please be (very) patient..."
|
||||||
|
unittest.main()
|
||||||
|
|
@ -13,7 +13,6 @@ class KernelTests(unittest.TestCase):
|
||||||
X = np.random.rand(5,5)
|
X = np.random.rand(5,5)
|
||||||
Y = np.ones((5,1))
|
Y = np.ones((5,1))
|
||||||
m = GPy.models.GP_regression(X,Y,K)
|
m = GPy.models.GP_regression(X,Y,K)
|
||||||
print m
|
|
||||||
self.assertTrue(m.checkgrad())
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
def test_coregionalisation(self):
|
def test_coregionalisation(self):
|
||||||
|
|
|
||||||
48
GPy/testing/sparse_gplvm_tests.py
Normal file
48
GPy/testing/sparse_gplvm_tests.py
Normal file
|
|
@ -0,0 +1,48 @@
|
||||||
|
# Copyright (c) 2012, Nicolo Fusi, James Hensman
|
||||||
|
# Licensed under the BSD 3-clause license (see LICENSE.txt)
|
||||||
|
|
||||||
|
import unittest
|
||||||
|
import numpy as np
|
||||||
|
import GPy
|
||||||
|
|
||||||
|
class sparse_GPLVMTests(unittest.TestCase):
|
||||||
|
def test_bias_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
k = GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.sparse_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
@unittest.skip('linear kernels do not have dKdiag_dX')
|
||||||
|
def test_linear_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
k = GPy.kern.linear(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.sparse_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
def test_rbf_kern(self):
|
||||||
|
N, M, Q, D = 10, 3, 2, 4
|
||||||
|
X = np.random.rand(N, Q)
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
K = k.K(X)
|
||||||
|
Y = np.random.multivariate_normal(np.zeros(N),K,D).T
|
||||||
|
k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001)
|
||||||
|
m = GPy.models.sparse_GPLVM(Y, Q, kernel = k, M=M)
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
m.randomize()
|
||||||
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
|
if __name__ == "__main__":
|
||||||
|
print "Running unit tests, please be (very) patient..."
|
||||||
|
unittest.main()
|
||||||
|
|
@ -192,17 +192,6 @@ class GradientTests(unittest.TestCase):
|
||||||
m.approximate_likelihood()
|
m.approximate_likelihood()
|
||||||
self.assertTrue(m.checkgrad())
|
self.assertTrue(m.checkgrad())
|
||||||
|
|
||||||
def test_warped_GP(self):
|
|
||||||
xmin, xmax = 1, 2.5*np.pi
|
|
||||||
b, C, SNR = 1, 0, 0.1
|
|
||||||
X = np.linspace(xmin, xmax, 500)
|
|
||||||
y = b*X + C + 1*np.sin(X)
|
|
||||||
y += 0.05*np.random.randn(len(X))
|
|
||||||
X, y = X[:, None], y[:, None]
|
|
||||||
m = GPy.models.warpedGP(X, y, warping_terms = 3)
|
|
||||||
m.constrain_positive('(tanh_a|tanh_b|rbf|white|bias)')
|
|
||||||
self.assertTrue(m.checkgrad())
|
|
||||||
|
|
||||||
|
|
||||||
if __name__ == "__main__":
|
if __name__ == "__main__":
|
||||||
print "Running unit tests, please be (very) patient..."
|
print "Running unit tests, please be (very) patient..."
|
||||||
|
|
|
||||||
|
|
@ -25,7 +25,7 @@ def fewerXticks(ax=None,divideby=2):
|
||||||
ax.set_xticks(ax.get_xticks()[::divideby])
|
ax.set_xticks(ax.get_xticks()[::divideby])
|
||||||
|
|
||||||
|
|
||||||
coloursHex = {\
|
colorsHex = {\
|
||||||
"Aluminium6":"#2e3436",\
|
"Aluminium6":"#2e3436",\
|
||||||
"Aluminium5":"#555753",\
|
"Aluminium5":"#555753",\
|
||||||
"Aluminium4":"#888a85",\
|
"Aluminium4":"#888a85",\
|
||||||
|
|
@ -54,9 +54,9 @@ coloursHex = {\
|
||||||
"mediumButter":"#edd400",\
|
"mediumButter":"#edd400",\
|
||||||
"darkButter":"#c4a000"}
|
"darkButter":"#c4a000"}
|
||||||
|
|
||||||
darkList = [coloursHex['darkBlue'],coloursHex['darkRed'],coloursHex['darkGreen'], coloursHex['darkOrange'], coloursHex['darkButter'], coloursHex['darkPurple'], coloursHex['darkChocolate'], coloursHex['Aluminium6']]
|
darkList = [colorsHex['darkBlue'],colorsHex['darkRed'],colorsHex['darkGreen'], colorsHex['darkOrange'], colorsHex['darkButter'], colorsHex['darkPurple'], colorsHex['darkChocolate'], colorsHex['Aluminium6']]
|
||||||
mediumList = [coloursHex['mediumBlue'], coloursHex['mediumRed'],coloursHex['mediumGreen'], coloursHex['mediumOrange'], coloursHex['mediumButter'], coloursHex['mediumPurple'], coloursHex['mediumChocolate'], coloursHex['Aluminium5']]
|
mediumList = [colorsHex['mediumBlue'], colorsHex['mediumRed'],colorsHex['mediumGreen'], colorsHex['mediumOrange'], colorsHex['mediumButter'], colorsHex['mediumPurple'], colorsHex['mediumChocolate'], colorsHex['Aluminium5']]
|
||||||
lightList = [coloursHex['lightBlue'], coloursHex['lightRed'],coloursHex['lightGreen'], coloursHex['lightOrange'], coloursHex['lightButter'], coloursHex['lightPurple'], coloursHex['lightChocolate'], coloursHex['Aluminium4']]
|
lightList = [colorsHex['lightBlue'], colorsHex['lightRed'],colorsHex['lightGreen'], colorsHex['lightOrange'], colorsHex['lightButter'], colorsHex['lightPurple'], colorsHex['lightChocolate'], colorsHex['Aluminium4']]
|
||||||
|
|
||||||
def currentDark():
|
def currentDark():
|
||||||
return darkList[-1]
|
return darkList[-1]
|
||||||
|
|
@ -76,85 +76,85 @@ def nextLight():
|
||||||
return lightList[-1]
|
return lightList[-1]
|
||||||
|
|
||||||
def reset():
|
def reset():
|
||||||
while not darkList[0]==coloursHex['darkBlue']:
|
while not darkList[0]==colorsHex['darkBlue']:
|
||||||
darkList.append(darkList.pop(0))
|
darkList.append(darkList.pop(0))
|
||||||
while not mediumList[0]==coloursHex['mediumBlue']:
|
while not mediumList[0]==colorsHex['mediumBlue']:
|
||||||
mediumList.append(mediumList.pop(0))
|
mediumList.append(mediumList.pop(0))
|
||||||
while not lightList[0]==coloursHex['lightBlue']:
|
while not lightList[0]==colorsHex['lightBlue']:
|
||||||
lightList.append(lightList.pop(0))
|
lightList.append(lightList.pop(0))
|
||||||
|
|
||||||
def setLightFigures():
|
def setLightFigures():
|
||||||
mpl.rcParams['axes.edgecolor']=coloursHex['Aluminium6']
|
mpl.rcParams['axes.edgecolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['axes.facecolor']=coloursHex['Aluminium2']
|
mpl.rcParams['axes.facecolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['axes.labelcolor']=coloursHex['Aluminium6']
|
mpl.rcParams['axes.labelcolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['figure.edgecolor']=coloursHex['Aluminium6']
|
mpl.rcParams['figure.edgecolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['figure.facecolor']=coloursHex['Aluminium2']
|
mpl.rcParams['figure.facecolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['grid.color']=coloursHex['Aluminium6']
|
mpl.rcParams['grid.color']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['savefig.edgecolor']=coloursHex['Aluminium2']
|
mpl.rcParams['savefig.edgecolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['savefig.facecolor']=coloursHex['Aluminium2']
|
mpl.rcParams['savefig.facecolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['text.color']=coloursHex['Aluminium6']
|
mpl.rcParams['text.color']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['xtick.color']=coloursHex['Aluminium6']
|
mpl.rcParams['xtick.color']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['ytick.color']=coloursHex['Aluminium6']
|
mpl.rcParams['ytick.color']=colorsHex['Aluminium6']
|
||||||
|
|
||||||
def setDarkFigures():
|
def setDarkFigures():
|
||||||
mpl.rcParams['axes.edgecolor']=coloursHex['Aluminium2']
|
mpl.rcParams['axes.edgecolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['axes.facecolor']=coloursHex['Aluminium6']
|
mpl.rcParams['axes.facecolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['axes.labelcolor']=coloursHex['Aluminium2']
|
mpl.rcParams['axes.labelcolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['figure.edgecolor']=coloursHex['Aluminium2']
|
mpl.rcParams['figure.edgecolor']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['figure.facecolor']=coloursHex['Aluminium6']
|
mpl.rcParams['figure.facecolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['grid.color']=coloursHex['Aluminium2']
|
mpl.rcParams['grid.color']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['savefig.edgecolor']=coloursHex['Aluminium6']
|
mpl.rcParams['savefig.edgecolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['savefig.facecolor']=coloursHex['Aluminium6']
|
mpl.rcParams['savefig.facecolor']=colorsHex['Aluminium6']
|
||||||
mpl.rcParams['text.color']=coloursHex['Aluminium2']
|
mpl.rcParams['text.color']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['xtick.color']=coloursHex['Aluminium2']
|
mpl.rcParams['xtick.color']=colorsHex['Aluminium2']
|
||||||
mpl.rcParams['ytick.color']=coloursHex['Aluminium2']
|
mpl.rcParams['ytick.color']=colorsHex['Aluminium2']
|
||||||
|
|
||||||
def hex2rgb(hexcolor):
|
def hex2rgb(hexcolor):
|
||||||
hexcolor = [hexcolor[1+2*i:1+2*(i+1)] for i in range(3)]
|
hexcolor = [hexcolor[1+2*i:1+2*(i+1)] for i in range(3)]
|
||||||
r,g,b = [int(n,16) for n in hexcolor]
|
r,g,b = [int(n,16) for n in hexcolor]
|
||||||
return (r,g,b)
|
return (r,g,b)
|
||||||
|
|
||||||
coloursRGB = dict([(k,hex2rgb(i)) for k,i in coloursHex.items()])
|
colorsRGB = dict([(k,hex2rgb(i)) for k,i in colorsHex.items()])
|
||||||
|
|
||||||
cdict_RB = {'red' :((0.,coloursRGB['mediumRed'][0]/256.,coloursRGB['mediumRed'][0]/256.),
|
cdict_RB = {'red' :((0.,colorsRGB['mediumRed'][0]/256.,colorsRGB['mediumRed'][0]/256.),
|
||||||
(.5,coloursRGB['mediumPurple'][0]/256.,coloursRGB['mediumPurple'][0]/256.),
|
(.5,colorsRGB['mediumPurple'][0]/256.,colorsRGB['mediumPurple'][0]/256.),
|
||||||
(1.,coloursRGB['mediumBlue'][0]/256.,coloursRGB['mediumBlue'][0]/256.)),
|
(1.,colorsRGB['mediumBlue'][0]/256.,colorsRGB['mediumBlue'][0]/256.)),
|
||||||
'green':((0.,coloursRGB['mediumRed'][1]/256.,coloursRGB['mediumRed'][1]/256.),
|
'green':((0.,colorsRGB['mediumRed'][1]/256.,colorsRGB['mediumRed'][1]/256.),
|
||||||
(.5,coloursRGB['mediumPurple'][1]/256.,coloursRGB['mediumPurple'][1]/256.),
|
(.5,colorsRGB['mediumPurple'][1]/256.,colorsRGB['mediumPurple'][1]/256.),
|
||||||
(1.,coloursRGB['mediumBlue'][1]/256.,coloursRGB['mediumBlue'][1]/256.)),
|
(1.,colorsRGB['mediumBlue'][1]/256.,colorsRGB['mediumBlue'][1]/256.)),
|
||||||
'blue':((0.,coloursRGB['mediumRed'][2]/256.,coloursRGB['mediumRed'][2]/256.),
|
'blue':((0.,colorsRGB['mediumRed'][2]/256.,colorsRGB['mediumRed'][2]/256.),
|
||||||
(.5,coloursRGB['mediumPurple'][2]/256.,coloursRGB['mediumPurple'][2]/256.),
|
(.5,colorsRGB['mediumPurple'][2]/256.,colorsRGB['mediumPurple'][2]/256.),
|
||||||
(1.,coloursRGB['mediumBlue'][2]/256.,coloursRGB['mediumBlue'][2]/256.))}
|
(1.,colorsRGB['mediumBlue'][2]/256.,colorsRGB['mediumBlue'][2]/256.))}
|
||||||
|
|
||||||
cdict_BGR = {'red' :((0.,coloursRGB['mediumBlue'][0]/256.,coloursRGB['mediumBlue'][0]/256.),
|
cdict_BGR = {'red' :((0.,colorsRGB['mediumBlue'][0]/256.,colorsRGB['mediumBlue'][0]/256.),
|
||||||
(.5,coloursRGB['mediumGreen'][0]/256.,coloursRGB['mediumGreen'][0]/256.),
|
(.5,colorsRGB['mediumGreen'][0]/256.,colorsRGB['mediumGreen'][0]/256.),
|
||||||
(1.,coloursRGB['mediumRed'][0]/256.,coloursRGB['mediumRed'][0]/256.)),
|
(1.,colorsRGB['mediumRed'][0]/256.,colorsRGB['mediumRed'][0]/256.)),
|
||||||
'green':((0.,coloursRGB['mediumBlue'][1]/256.,coloursRGB['mediumBlue'][1]/256.),
|
'green':((0.,colorsRGB['mediumBlue'][1]/256.,colorsRGB['mediumBlue'][1]/256.),
|
||||||
(.5,coloursRGB['mediumGreen'][1]/256.,coloursRGB['mediumGreen'][1]/256.),
|
(.5,colorsRGB['mediumGreen'][1]/256.,colorsRGB['mediumGreen'][1]/256.),
|
||||||
(1.,coloursRGB['mediumRed'][1]/256.,coloursRGB['mediumRed'][1]/256.)),
|
(1.,colorsRGB['mediumRed'][1]/256.,colorsRGB['mediumRed'][1]/256.)),
|
||||||
'blue':((0.,coloursRGB['mediumBlue'][2]/256.,coloursRGB['mediumBlue'][2]/256.),
|
'blue':((0.,colorsRGB['mediumBlue'][2]/256.,colorsRGB['mediumBlue'][2]/256.),
|
||||||
(.5,coloursRGB['mediumGreen'][2]/256.,coloursRGB['mediumGreen'][2]/256.),
|
(.5,colorsRGB['mediumGreen'][2]/256.,colorsRGB['mediumGreen'][2]/256.),
|
||||||
(1.,coloursRGB['mediumRed'][2]/256.,coloursRGB['mediumRed'][2]/256.))}
|
(1.,colorsRGB['mediumRed'][2]/256.,colorsRGB['mediumRed'][2]/256.))}
|
||||||
|
|
||||||
|
|
||||||
cdict_Alu = {'red' :((0./5,coloursRGB['Aluminium1'][0]/256.,coloursRGB['Aluminium1'][0]/256.),
|
cdict_Alu = {'red' :((0./5,colorsRGB['Aluminium1'][0]/256.,colorsRGB['Aluminium1'][0]/256.),
|
||||||
(1./5,coloursRGB['Aluminium2'][0]/256.,coloursRGB['Aluminium2'][0]/256.),
|
(1./5,colorsRGB['Aluminium2'][0]/256.,colorsRGB['Aluminium2'][0]/256.),
|
||||||
(2./5,coloursRGB['Aluminium3'][0]/256.,coloursRGB['Aluminium3'][0]/256.),
|
(2./5,colorsRGB['Aluminium3'][0]/256.,colorsRGB['Aluminium3'][0]/256.),
|
||||||
(3./5,coloursRGB['Aluminium4'][0]/256.,coloursRGB['Aluminium4'][0]/256.),
|
(3./5,colorsRGB['Aluminium4'][0]/256.,colorsRGB['Aluminium4'][0]/256.),
|
||||||
(4./5,coloursRGB['Aluminium5'][0]/256.,coloursRGB['Aluminium5'][0]/256.),
|
(4./5,colorsRGB['Aluminium5'][0]/256.,colorsRGB['Aluminium5'][0]/256.),
|
||||||
(5./5,coloursRGB['Aluminium6'][0]/256.,coloursRGB['Aluminium6'][0]/256.)),
|
(5./5,colorsRGB['Aluminium6'][0]/256.,colorsRGB['Aluminium6'][0]/256.)),
|
||||||
'green' :((0./5,coloursRGB['Aluminium1'][1]/256.,coloursRGB['Aluminium1'][1]/256.),
|
'green' :((0./5,colorsRGB['Aluminium1'][1]/256.,colorsRGB['Aluminium1'][1]/256.),
|
||||||
(1./5,coloursRGB['Aluminium2'][1]/256.,coloursRGB['Aluminium2'][1]/256.),
|
(1./5,colorsRGB['Aluminium2'][1]/256.,colorsRGB['Aluminium2'][1]/256.),
|
||||||
(2./5,coloursRGB['Aluminium3'][1]/256.,coloursRGB['Aluminium3'][1]/256.),
|
(2./5,colorsRGB['Aluminium3'][1]/256.,colorsRGB['Aluminium3'][1]/256.),
|
||||||
(3./5,coloursRGB['Aluminium4'][1]/256.,coloursRGB['Aluminium4'][1]/256.),
|
(3./5,colorsRGB['Aluminium4'][1]/256.,colorsRGB['Aluminium4'][1]/256.),
|
||||||
(4./5,coloursRGB['Aluminium5'][1]/256.,coloursRGB['Aluminium5'][1]/256.),
|
(4./5,colorsRGB['Aluminium5'][1]/256.,colorsRGB['Aluminium5'][1]/256.),
|
||||||
(5./5,coloursRGB['Aluminium6'][1]/256.,coloursRGB['Aluminium6'][1]/256.)),
|
(5./5,colorsRGB['Aluminium6'][1]/256.,colorsRGB['Aluminium6'][1]/256.)),
|
||||||
'blue' :((0./5,coloursRGB['Aluminium1'][2]/256.,coloursRGB['Aluminium1'][2]/256.),
|
'blue' :((0./5,colorsRGB['Aluminium1'][2]/256.,colorsRGB['Aluminium1'][2]/256.),
|
||||||
(1./5,coloursRGB['Aluminium2'][2]/256.,coloursRGB['Aluminium2'][2]/256.),
|
(1./5,colorsRGB['Aluminium2'][2]/256.,colorsRGB['Aluminium2'][2]/256.),
|
||||||
(2./5,coloursRGB['Aluminium3'][2]/256.,coloursRGB['Aluminium3'][2]/256.),
|
(2./5,colorsRGB['Aluminium3'][2]/256.,colorsRGB['Aluminium3'][2]/256.),
|
||||||
(3./5,coloursRGB['Aluminium4'][2]/256.,coloursRGB['Aluminium4'][2]/256.),
|
(3./5,colorsRGB['Aluminium4'][2]/256.,colorsRGB['Aluminium4'][2]/256.),
|
||||||
(4./5,coloursRGB['Aluminium5'][2]/256.,coloursRGB['Aluminium5'][2]/256.),
|
(4./5,colorsRGB['Aluminium5'][2]/256.,colorsRGB['Aluminium5'][2]/256.),
|
||||||
(5./5,coloursRGB['Aluminium6'][2]/256.,coloursRGB['Aluminium6'][2]/256.))}
|
(5./5,colorsRGB['Aluminium6'][2]/256.,colorsRGB['Aluminium6'][2]/256.))}
|
||||||
# cmap_Alu = mpl.colors.LinearSegmentedColormap('TangoAluminium',cdict_Alu,256)
|
# cmap_Alu = mpl.colors.LinearSegmentedColormap('TangoAluminium',cdict_Alu,256)
|
||||||
# cmap_BGR = mpl.colors.LinearSegmentedColormap('TangoRedBlue',cdict_BGR,256)
|
# cmap_BGR = mpl.colors.LinearSegmentedColormap('TangoRedBlue',cdict_BGR,256)
|
||||||
# cmap_RB = mpl.colors.LinearSegmentedColormap('TangoRedBlue',cdict_RB,256)
|
# cmap_RB = mpl.colors.LinearSegmentedColormap('TangoRedBlue',cdict_RB,256)
|
||||||
|
|
|
||||||
|
|
@ -46,7 +46,7 @@ def oil_100(seed=default_seed):
|
||||||
return {'X': X, 'Y': Y, 'info': "Subsample of the oil data extracting 100 values randomly without replacement."}
|
return {'X': X, 'Y': Y, 'info': "Subsample of the oil data extracting 100 values randomly without replacement."}
|
||||||
|
|
||||||
def pumadyn(seed=default_seed):
|
def pumadyn(seed=default_seed):
|
||||||
# Data is variance 1, no need to normalise.
|
# Data is variance 1, no need to normalize.
|
||||||
data = np.loadtxt(os.path.join(data_path, 'pumadyn-32nm/Dataset.data.gz'))
|
data = np.loadtxt(os.path.join(data_path, 'pumadyn-32nm/Dataset.data.gz'))
|
||||||
indices = np.random.permutation(data.shape[0])
|
indices = np.random.permutation(data.shape[0])
|
||||||
indicesTrain = indices[0:7168]
|
indicesTrain = indices[0:7168]
|
||||||
|
|
|
||||||
|
|
@ -11,7 +11,7 @@ import re
|
||||||
import pdb
|
import pdb
|
||||||
import cPickle
|
import cPickle
|
||||||
import types
|
import types
|
||||||
import scipy.lib.lapack.flapack
|
#import scipy.lib.lapack.flapack
|
||||||
import scipy as sp
|
import scipy as sp
|
||||||
|
|
||||||
def mdot(*args):
|
def mdot(*args):
|
||||||
|
|
@ -101,7 +101,7 @@ def chol_inv(L):
|
||||||
|
|
||||||
"""
|
"""
|
||||||
|
|
||||||
return linalg.flapack.dtrtri(L, lower = True)[0]
|
return linalg.lapack.flapack.dtrtri(L, lower = True)[0]
|
||||||
|
|
||||||
|
|
||||||
def multiple_pdinv(A):
|
def multiple_pdinv(A):
|
||||||
|
|
@ -118,7 +118,7 @@ def multiple_pdinv(A):
|
||||||
N = A.shape[-1]
|
N = A.shape[-1]
|
||||||
chols = [jitchol(A[:,:,i]) for i in range(N)]
|
chols = [jitchol(A[:,:,i]) for i in range(N)]
|
||||||
halflogdets = [np.sum(np.log(np.diag(L[0]))) for L in chols]
|
halflogdets = [np.sum(np.log(np.diag(L[0]))) for L in chols]
|
||||||
invs = [linalg.flapack.dpotri(L[0],True)[0] for L in chols]
|
invs = [linalg.lapack.flapack.dpotri(L[0],True)[0] for L in chols]
|
||||||
invs = [np.triu(I)+np.triu(I,1).T for I in invs]
|
invs = [np.triu(I)+np.triu(I,1).T for I in invs]
|
||||||
return np.dstack(invs),np.array(halflogdets)
|
return np.dstack(invs),np.array(halflogdets)
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -6,7 +6,7 @@ import Tango
|
||||||
import pylab as pb
|
import pylab as pb
|
||||||
import numpy as np
|
import numpy as np
|
||||||
|
|
||||||
def gpplot(x,mu,lower,upper,edgecol=Tango.coloursHex['darkBlue'],fillcol=Tango.coloursHex['lightBlue'],axes=None,**kwargs):
|
def gpplot(x,mu,lower,upper,edgecol=Tango.colorsHex['darkBlue'],fillcol=Tango.colorsHex['lightBlue'],axes=None,**kwargs):
|
||||||
if axes is None:
|
if axes is None:
|
||||||
axes = pb.gca()
|
axes = pb.gca()
|
||||||
mu = mu.flatten()
|
mu = mu.flatten()
|
||||||
|
|
|
||||||
BIN
doc/Figures/tick.png
Normal file
BIN
doc/Figures/tick.png
Normal file
Binary file not shown.
|
After Width: | Height: | Size: 175 B |
|
|
@ -73,6 +73,22 @@ examples Package
|
||||||
:undoc-members:
|
:undoc-members:
|
||||||
:show-inheritance:
|
:show-inheritance:
|
||||||
|
|
||||||
|
:mod:`tuto_GP_regression` Module
|
||||||
|
--------------------------------
|
||||||
|
|
||||||
|
.. automodule:: GPy.examples.tuto_GP_regression
|
||||||
|
:members:
|
||||||
|
:undoc-members:
|
||||||
|
:show-inheritance:
|
||||||
|
|
||||||
|
:mod:`tuto_kernel_overview` Module
|
||||||
|
----------------------------------
|
||||||
|
|
||||||
|
.. automodule:: GPy.examples.tuto_kernel_overview
|
||||||
|
:members:
|
||||||
|
:undoc-members:
|
||||||
|
:show-inheritance:
|
||||||
|
|
||||||
:mod:`uncertain_input_GP_regression_demo` Module
|
:mod:`uncertain_input_GP_regression_demo` Module
|
||||||
------------------------------------------------
|
------------------------------------------------
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -49,6 +49,14 @@ kern Package
|
||||||
:undoc-members:
|
:undoc-members:
|
||||||
:show-inheritance:
|
:show-inheritance:
|
||||||
|
|
||||||
|
:mod:`coregionalise` Module
|
||||||
|
---------------------------
|
||||||
|
|
||||||
|
.. automodule:: GPy.kern.coregionalise
|
||||||
|
:members:
|
||||||
|
:undoc-members:
|
||||||
|
:show-inheritance:
|
||||||
|
|
||||||
:mod:`exponential` Module
|
:mod:`exponential` Module
|
||||||
-------------------------
|
-------------------------
|
||||||
|
|
||||||
|
|
@ -113,18 +121,18 @@ kern Package
|
||||||
:undoc-members:
|
:undoc-members:
|
||||||
:show-inheritance:
|
:show-inheritance:
|
||||||
|
|
||||||
:mod:`product` Module
|
:mod:`prod` Module
|
||||||
---------------------
|
------------------
|
||||||
|
|
||||||
.. automodule:: GPy.kern.product
|
.. automodule:: GPy.kern.prod
|
||||||
:members:
|
:members:
|
||||||
:undoc-members:
|
:undoc-members:
|
||||||
:show-inheritance:
|
:show-inheritance:
|
||||||
|
|
||||||
:mod:`product_orthogonal` Module
|
:mod:`prod_orthogonal` Module
|
||||||
--------------------------------
|
-----------------------------
|
||||||
|
|
||||||
.. automodule:: GPy.kern.product_orthogonal
|
.. automodule:: GPy.kern.prod_orthogonal
|
||||||
:members:
|
:members:
|
||||||
:undoc-members:
|
:undoc-members:
|
||||||
:show-inheritance:
|
:show-inheritance:
|
||||||
|
|
@ -145,6 +153,14 @@ kern Package
|
||||||
:undoc-members:
|
:undoc-members:
|
||||||
:show-inheritance:
|
:show-inheritance:
|
||||||
|
|
||||||
|
:mod:`symmetric` Module
|
||||||
|
-----------------------
|
||||||
|
|
||||||
|
.. automodule:: GPy.kern.symmetric
|
||||||
|
:members:
|
||||||
|
:undoc-members:
|
||||||
|
:show-inheritance:
|
||||||
|
|
||||||
:mod:`sympykern` Module
|
:mod:`sympykern` Module
|
||||||
-----------------------
|
-----------------------
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -3,15 +3,45 @@
|
||||||
List of implemented kernels
|
List of implemented kernels
|
||||||
***************************
|
***************************
|
||||||
|
|
||||||
The :math:`\checkmark` symbol represents the functions that have been implemented for each kernel.
|
The following table shows the implemented kernels in GPy and gives the details of the implemented function for each kernel.
|
||||||
|
|
||||||
.. |tick|
|
==================== =========== ====== ======= =========== =============== ======= =========== ====== ====== =======
|
||||||
|
NAME get/set K Kdiag dK_dtheta dKdiag_dtheta dK_dX dKdiag_dX psi0 psi1 psi2
|
||||||
|
==================== =========== ====== ======= =========== =============== ======= =========== ====== ====== =======
|
||||||
|
bias |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
Brownian |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
exponential |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
finite_dimensional |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
linear |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
Matern32 |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
Matern52 |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
periodic_exponential |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
periodic_Matern32 |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
periodic_Matern52 |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
rbf |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
spline |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
-------------------- ----------- ------ ------- ----------- --------------- ------- ----------- ------ ------ -------
|
||||||
|
white |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick| |tick|
|
||||||
|
==================== =========== ====== ======= =========== =============== ======= =========== ====== ====== =======
|
||||||
|
|
||||||
.. |tick| image:: tick.png
|
Depending on the use, all functions may not be required
|
||||||
|
|
||||||
|
* ``get/set, K, Kdiag``: compulsory
|
||||||
|
* ``dK_dtheta``: necessary to optimize the model
|
||||||
|
* ``dKdiag_dtheta``: sparse models, BGPLVM, GPs with uncertain inputs
|
||||||
|
* ``dK_dX``: sparse models, GPLVM, BGPLVM, GPs with uncertain inputs
|
||||||
|
* ``dKdiag_dX``: sparse models, BGPLVM, GPs with uncertain inputs
|
||||||
|
* ``psi0, psi1, psi2``: BGPLVM, GPs with uncertain inputs
|
||||||
|
|
||||||
====== =========== === ======= =========== =============== ======= =========== ====== ====== =======
|
.. |tick| image:: Figures/tick.png
|
||||||
NAME get/set K Kdiag dK_dtheta dKdiag_dtheta dK_dX dKdiag_dX psi0 psi1 psi2
|
|
||||||
====== =========== === ======= =========== =============== ======= =========== ====== ====== =======
|
|
||||||
rbf \\checkmark y
|
|
||||||
====== =========== === ======= =========== =============== ======= =========== ====== ====== =======
|
|
||||||
|
|
|
||||||
|
|
@ -2,7 +2,7 @@
|
||||||
Gaussian process regression tutorial
|
Gaussian process regression tutorial
|
||||||
*************************************
|
*************************************
|
||||||
|
|
||||||
We will see in this tutorial the basics for building a 1 dimensional and a 2 dimensional Gaussian process regression model, also known as a kriging model. The code shown in this tutorial can be found without the comments at GPy/examples/tuto_GP_regression.py.
|
We will see in this tutorial the basics for building a 1 dimensional and a 2 dimensional Gaussian process regression model, also known as a kriging model. The code shown in this tutorial can be obtained at GPy/examples/tutorials.py, or by running ``GPy.examples.tutorials.tuto_GP_regression()``.
|
||||||
|
|
||||||
We first import the libraries we will need: ::
|
We first import the libraries we will need: ::
|
||||||
|
|
||||||
|
|
|
||||||
60
doc/tuto_interacting_with_models.rst
Normal file
60
doc/tuto_interacting_with_models.rst
Normal file
|
|
@ -0,0 +1,60 @@
|
||||||
|
*************************************
|
||||||
|
Interacting with models
|
||||||
|
*************************************
|
||||||
|
|
||||||
|
The GPy model class has a set of features which are designed to make it simple to explore the parameter space of the model. By default, the scipy optimisers are used to fit GPy models (via model.optimize()), for which we provide mechanisms for 'free' optimisation: GPy can ensure that naturally positive parameters (such as variances) remain positive. But these mechanisms are much more powerful than simple reparameterisation, as we shall see.
|
||||||
|
|
||||||
|
All of the examples included in GPy return an instance of a model class. We'll use GPy.examples.?? as an example::
|
||||||
|
|
||||||
|
import pylab as pb
|
||||||
|
pb.ion()
|
||||||
|
import GPy
|
||||||
|
m = GPy.examples.??
|
||||||
|
|
||||||
|
Examining the model using print
|
||||||
|
===============================
|
||||||
|
To see the current state of the model parameters, and the model's (marginal) likelihood just print the model::
|
||||||
|
print m
|
||||||
|
|
||||||
|
?? output
|
||||||
|
|
||||||
|
Getting the model's likelihood and gradients
|
||||||
|
===========================================
|
||||||
|
foobar
|
||||||
|
|
||||||
|
Setting and fetching parameters by name
|
||||||
|
=======================================
|
||||||
|
foobar
|
||||||
|
|
||||||
|
Constraining and optimising the model
|
||||||
|
=====================================
|
||||||
|
A simple task in GPy is to ensure that the models' variances remain positive during optimisation. the models class has a function called constrain_positive(), which accepts a regex string as above. To constrain the models' variance to be positive::
|
||||||
|
m.constrain_positive('variance')
|
||||||
|
print m
|
||||||
|
|
||||||
|
Now we see that the variance of the model is constrained to be postive. GPy handles the effective change of gradients: see how m.objective_gradients has changed approriately
|
||||||
|
|
||||||
|
|
||||||
|
For convenience, we also provide a catch all function which ensures that anything which appears to require positivity is constrianed appropriately::
|
||||||
|
m.ensure_default_constraints()
|
||||||
|
|
||||||
|
|
||||||
|
Fixing parameters
|
||||||
|
=================
|
||||||
|
|
||||||
|
|
||||||
|
Tying Parameters
|
||||||
|
================
|
||||||
|
|
||||||
|
Bounding parameters
|
||||||
|
===================
|
||||||
|
|
||||||
|
|
||||||
|
Further Reading
|
||||||
|
===============
|
||||||
|
All of the mechansiams for dealing with parameters are baked right into GPy.core.model, from which all of the classes in GPy.models inherrit. To learn how to construct your own model, you might want to read ??link?? creating_new_models.
|
||||||
|
|
||||||
|
By deafult, GPy uses the tnc optimizer (from scipy.optimize.tnc). To use other optimisers, and to control the setting of those optimisers, as well as other funky features like automated restarts and diagnostics, you can read the optimization tutorial ??link??.
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
@ -2,7 +2,7 @@
|
||||||
****************************
|
****************************
|
||||||
tutorial : A kernel overview
|
tutorial : A kernel overview
|
||||||
****************************
|
****************************
|
||||||
The aim of this tutorial is to give a better understanding of the kernel objects in GPy and to list the ones that are already implemented. The code shown in this tutorial can be found without the comments at GPy/examples/tuto_kernel_overview.py.
|
The aim of this tutorial is to give a better understanding of the kernel objects in GPy and to list the ones that are already implemented. The code shown in this tutorial can be obtained at GPy/examples/tutorials.py or by running ``GPy.examples.tutorials.tuto_kernel_overview()``.
|
||||||
|
|
||||||
First we import the libraries we will need ::
|
First we import the libraries we will need ::
|
||||||
|
|
||||||
|
|
@ -39,7 +39,7 @@ return::
|
||||||
Implemented kernels
|
Implemented kernels
|
||||||
===================
|
===================
|
||||||
|
|
||||||
Many kernels are already implemented in GPy. A comprehensive list can be found `here <kernel_implementation.html>`_ . The following figure gives a summary of most of them:
|
Many kernels are already implemented in GPy. A comprehensive list can be found `here <kernel_implementation.html>`_ and the following figure gives a summary of most of them:
|
||||||
|
|
||||||
.. figure:: Figures/tuto_kern_overview_allkern.png
|
.. figure:: Figures/tuto_kern_overview_allkern.png
|
||||||
:align: center
|
:align: center
|
||||||
|
|
|
||||||
8
setup.py
8
setup.py
|
|
@ -3,8 +3,6 @@
|
||||||
|
|
||||||
import os
|
import os
|
||||||
from setuptools import setup
|
from setuptools import setup
|
||||||
#from numpy.distutils.core import Extension, setup
|
|
||||||
#from sphinx.setup_command import BuildDoc
|
|
||||||
|
|
||||||
# Version number
|
# Version number
|
||||||
version = '0.1.3'
|
version = '0.1.3'
|
||||||
|
|
@ -14,12 +12,12 @@ def read(fname):
|
||||||
|
|
||||||
setup(name = 'GPy',
|
setup(name = 'GPy',
|
||||||
version = version,
|
version = version,
|
||||||
author = 'James Hensman, Nicolo Fusi, Ricardo Andrade, Nicolas Durrande, Alan Saul, Neil D. Lawrence',
|
author = read('AUTHORS.txt'),
|
||||||
author_email = "james.hensman@gmail.com",
|
author_email = "james.hensman@gmail.com",
|
||||||
description = ("The Gaussian Process Toolbox"),
|
description = ("The Gaussian Process Toolbox"),
|
||||||
license = "BSD 3-clause",
|
license = "BSD 3-clause",
|
||||||
keywords = "machine-learning gaussian-processes kernels",
|
keywords = "machine-learning gaussian-processes kernels",
|
||||||
url = "http://ml.sheffield.ac.uk/GPy/",
|
url = "http://sheffieldml.github.com/GPy/",
|
||||||
packages = ['GPy', 'GPy.core', 'GPy.kern', 'GPy.util', 'GPy.models', 'GPy.inference', 'GPy.examples', 'GPy.likelihoods'],
|
packages = ['GPy', 'GPy.core', 'GPy.kern', 'GPy.util', 'GPy.models', 'GPy.inference', 'GPy.examples', 'GPy.likelihoods'],
|
||||||
package_dir={'GPy': 'GPy'},
|
package_dir={'GPy': 'GPy'},
|
||||||
package_data = {'GPy': ['GPy/examples']},
|
package_data = {'GPy': ['GPy/examples']},
|
||||||
|
|
@ -34,7 +32,5 @@ setup(name = 'GPy',
|
||||||
#setup_requires=['sphinx'],
|
#setup_requires=['sphinx'],
|
||||||
#cmdclass = {'build_sphinx': BuildDoc},
|
#cmdclass = {'build_sphinx': BuildDoc},
|
||||||
classifiers=[
|
classifiers=[
|
||||||
"Development Status :: 1 - Alpha",
|
|
||||||
"Topic :: Machine Learning",
|
|
||||||
"License :: OSI Approved :: BSD License"],
|
"License :: OSI Approved :: BSD License"],
|
||||||
)
|
)
|
||||||
|
|
|
||||||
Loading…
Add table
Add a link
Reference in a new issue