apunkt/GPy
1
0
Fork 0
mirror of https://github.com/SheffieldML/GPy.git synced 2026-06-11 15:15:15 +02:00
Code Issues Projects Releases Packages Wiki Activity Actions

weaved coregionalise. much performance gained

Browse source
This commit is contained in:
James Hensman 2013-04-28 21:37:36 +01:00
parent 52ba8e4ba3
commit 7d9352c733
4 changed files with 70 additions and 11 deletions
Download patch file Download diff file Expand all files Collapse all files

62
GPy/kern/coregionalise.py
View file

@ -5,10 +5,11 @@ from kernpart import kernpart
import numpy as np
from GPy.util.linalg import mdot, pdinv
import pdb
from scipy import weave
class coregionalise(kernpart):
"""
Kernel for Intrisec Corregionalization Models
Kernel for Intrinsic Corregionalization Models
"""
def __init__(self,Nout,R=1, W=None, kappa=None):
self.D = 1
@ -42,19 +43,70 @@ class coregionalise(kernpart):
def K(self,index,index2,target):
index = np.asarray(index,dtype=np.int)
#here's the old code (numpy)
#if index2 is None:
#index2 = index
#else:
#index2 = np.asarray(index2,dtype=np.int)
#false_target = target.copy()
#ii,jj = np.meshgrid(index,index2)
#ii,jj = ii.T, jj.T
#false_target += self.B[ii,jj]
if index2 is None:
index2 = index
code="""
for(int i=0;i<N; i++){
target[i+i*N] += B[index[i]+Nout*index[i]];
for(int j=0; j<i; j++){
target[j+i*N] += B[index[i]+Nout*index[j]];
target[i+j*N] += target[j+i*N];
}
}
"""
N,B,Nout = index.size, self.B, self.Nout
weave.inline(code,['target','index','N','B','Nout'])
else:
index2 = np.asarray(index2,dtype=np.int)
ii,jj = np.meshgrid(index,index2)
ii,jj = ii.T, jj.T
target += self.B[ii,jj]
code="""
for(int i=0;i<M; i++){
for(int j=0; j<N; j++){
target[i+j*M] += B[Nout*index[j]+index2[i]];
}
}
"""
N,M,B,Nout = index.size,index2.size, self.B, self.Nout
weave.inline(code,['target','index','index2','N','M','B','Nout'])
def Kdiag(self,index,target):
target += np.diag(self.B)[np.asarray(index,dtype=np.int).flatten()]
def dK_dtheta(self,dL_dK,index,index2,target):
index = np.asarray(index,dtype=np.int)
dL_dK_small = np.zeros_like(self.B)
if index2 is None:
index2 = index
else:
index2 = np.asarray(index2,dtype=np.int)
code="""
for(int i=0; i<M; i++){
for(int j=0; j<N; j++){
dL_dK_small[index[j] + Nout*index2[i]] += dL_dK[i+j*M];
}
}
"""
N, M, Nout = index.size, index2.size, self.Nout
weave.inline(code, ['N','M','Nout','dL_dK','dL_dK_small','index','index2'])
dkappa = np.diag(dL_dK_small)
dL_dK_small += dL_dK_small.T
dW = (self.W[:,None,:]*dL_dK_small[:,:,None]).sum(0)
target += np.hstack([dW.flatten(),dkappa])
def dK_dtheta_old(self,dL_dK,index,index2,target):
if index2 is None:
index2 = index
else:

1
GPy/kern/kern.py
View file

@ -9,6 +9,7 @@ from kernpart import kernpart
import itertools
from prod_orthogonal import prod_orthogonal
from prod import prod
from ..util.linalg import symmetrify
class kern(parameterised):
def __init__(self, D, parts=[], input_slices=None):

9
GPy/kern/prod.py
View file

@ -40,9 +40,12 @@ class prod(kernpart):
def K(self,X,X2,target):
"""Compute the covariance matrix between X and X2."""
if X2 is None: X2 = X
target1 = np.zeros((X.shape[0],X2.shape[0]))
target2 = np.zeros((X.shape[0],X2.shape[0]))
if X2 is None:
target1 = np.zeros((X.shape[0],X2.shape[0]))
target2 = np.zeros((X.shape[0],X2.shape[0]))
else:
target1 = np.zeros((X.shape[0],X.shape[0]))
target2 = np.zeros((X.shape[0],X.shape[0]))
self.k1.K(X,X2,target1)
self.k2.K(X,X2,target2)
target += target1 * target2

9
GPy/kern/prod_orthogonal.py
View file

@ -39,11 +39,14 @@ class prod_orthogonal(kernpart):
def K(self,X,X2,target):
"""Compute the covariance matrix between X and X2."""
if X2 is None: X2 = X
target1 = np.zeros_like(target)
target2 = np.zeros_like(target)
self.k1.K(X[:,:self.k1.D],X2[:,:self.k1.D],target1)
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],target2)
if X2 is None:
self.k1.K(X[:,:self.k1.D],None,target1)
self.k2.K(X[:,self.k1.D:],None,target2)
else:
self.k1.K(X[:,:self.k1.D],X2[:,:self.k1.D],target1)
self.k2.K(X[:,self.k1.D:],X2[:,self.k1.D:],target2)
target += target1 * target2
def dK_dtheta(self,dL_dK,X,X2,target):