<commit_before># Gaussian process with derivatives in c
#LK_only can only be used for log likelihood or posterior log likelihood
#GPcore takes a single hyperparameter value or a set. With a set the infer_ methods produce individual inferences while infer_post methods return the posterior
#cython: profile=True
#

__author__ = "mark"
__date__ = "$22-Nov-2015 21:27:19$"

#TODO draw doesn't work correclty on a GP created with no data
import numpy as np
cimport numpy as np
import scipy as sp
import ctypes as ct
import os
import sys
from copy import deepcopy as dc
from scipy.stats import norm as norms
from libc.math cimport log10, log, isnan
#print os.path.join(os.path.split(__file__)[0],'../../dist/Release/GNU-Linux/libGPshared.so')
from . import __file__ as fl


libGP = ct.cdll.LoadLibrary(os.path.join(os.path.split(fl)[0],'../cproj/libcproj.so')) #path to c-shared library
#print libGP
libGP.k.restype = ct.c_double

ctpd = ct.POINTER(ct.c_double)
cint = ct.c_int


class GP_LKonly:
    def __init__(self, X_s, Y_s, S_s, D_s, kf):
        cdef int n, D, i
        n ,D = X_s.shape
        self.hyp = kf.hyp
        R = ct.c_double()
        Dc = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_s]
        libGP.newGP_LKonly(cint(D),cint(n),X_s.ctypes.data_as(ctpd),Y_s.ctypes.data_as(ctpd),S_s.ctypes.data_as(ctpd),(cint*len(Dc))(*Dc), cint(kf.Kindex), kf.hyp.ctypes.data_as(ctpd),ct.byref(R))
        self.l = R.value
        
        return
    
    def llk(self):
        #log likelihood of data given hyperparametrs
        return self.l
    
    def plk(self,lm,ls):
        #log likelihood given lognormal prior over hyperparameters
        tmp = 0.

        for i,h in enumerate(self.hyp):
            tmp -= 0.5*((log10(h)-lm[i])**2)/ls[i]**2
        
        return self.l+tmp

class GPcore:
    def __init__(self, X_s, Y_s, S_s, D_s, kf):
        #print [X_s, Y_s, S_s, D_s, kf]
        if type(kf) is kernel:
            self.size = 1
            kf = [kf]
        else:
            self.size = len(kf)
        allhyp = sp.hstack([k.hyp for k in kf])
        self.kf=kf
        [self.n ,self.D] = X_s.shape
        Dx = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_s]
        self.s = libGP.newGP_hypset(cint(self.D),cint(self.n),cint(kf[0].Kindex),X_s.ctypes.data_as(ctpd),Y_s.ctypes.data_as(ctpd),S_s.ctypes.data_as(ctpd),(cint*len(Dx))(*Dx),allhyp.ctypes.data_as(ctpd),cint(self.size))
        self.Y_s=Y_s
        
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_s]
        #print self.get_cho()
        libGP.presolv(self.s,cint(self.size))
        
        return
    
    def __del__(self):
        #print "here"
        libGP.killGP(self.s,cint(self.size))
        #print "also here"
        return
    
    def printc(self):
        print self.size
        libGP.ping(self.s, cint(self.size))
        return
    def get_cho(self):
        C=sp.empty([self.n,self.n*self.size])
        libGP.get_cho(self.s,cint(self.size),C.ctypes.data_as(ctpd))
        return C
    def infer_m(self,X_i,D_i):
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_i]
        R=sp.vstack([sp.empty(ns)]*self.size)
        libGP.infer_m(self.s, cint(self.size), ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
        return R
    
    def infer_m_partial(self,X_i,D_i,ki,hyp):
        
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_i]
        R=sp.vstack([sp.empty(ns)]*1)
        
        #libGP.infer_m_partial(self.s, cint(ki),hyp.ctypes.data_as(ctpd),ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
        libGP.infer_m_partial(self.s,cint(ki),hyp.ctypes.data_as(ctpd),ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
            
        
        return R
    
    def infer_m_post(self,X_i,D_i):
        ns=X_i.shape[0]
        R = self.infer_m(X_i,D_i)
        
        return sp.mean(R,axis=0).reshape([1,ns])
    
    
    def infer_full(self,X_i,D_i):
        cdef int ns,j,i
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**j for j in x])) for x in D_i]
        R=sp.vstack([sp.empty([ns+1,ns])]*self.size)
        libGP.infer_full(self.s, cint(self.size), ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
        m = sp.vstack([R[i*(ns+1),:] for i in range(self.size)])
        V = sp.vstack([R[(ns+1)*i+1:(ns+1)*(i+1),:] for i in range(self.size)])
        return [m,V]
    
    def infer_full_post(self,X_i,D_i):
        class MJMError(Exception):
            pass
        [m,V] = self.infer_full(X_i,D_i)
        cdef int i,ns
        ns=X_i.shape[0]
        cv = sp.zeros([ns,ns])

        for i in range(self.size):
            cv+=V[ns*i:ns*(i+1),:]
        cv= cv/self.size + sp.cov(m,rowvar=0,bias=1)
        return [sp.mean(m,axis=0).reshape([1,ns]),cv]
    
    def infer_diag(self,X_i,D_i):
        cdef int i,j,ns
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**j for j in x])) for x in D_i]
        R=sp.vstack([sp.empty([2,ns])]*self.size)
        libGP.infer_diag(self.s,cint(self.size), ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))

        m = sp.vstack([R[i*2,:] for i in range(self.size)])
        V = sp.vstack([R[i*2+1,:] for i in range(self.size)])
        return [m,V]
    
    def infer_diag_post(self,X_ii,D_i):
        cdef int ns
        X_i = dc(X_ii)
        ns = len(D_i)
        
        X_i.resize([ns,self.D])
        [m,V] = self.infer_diag(X_i,D_i)
        if sp.amin(V)<=-0.:
            class MJMError(Exception):
                pass
            print "negative/eq variance"
            print [m,V,X_i,D_i]
            print "_______________"
            #self.printc()
            raise(MJMError)
        if sp.amin(sp.var(m,axis=0))<-0.:
            class MJMError(Exception):
                pass
            print "negativevar of mean"
            print X_i.shape
            print [m,V,sp.var(m,axis=0),X_i,D_i]
            print "_______________"
            #self.printc()
            raise(MJMError)
        
        return [sp.mean(m,axis=0).reshape([1,ns]),(sp.mean(V,axis=0)+sp.var(m,axis=0)).reshape([1,ns])]
        
    
    def draw(self,X_i,D_i,m):
        #make m draws at X_i Nd, X, D, R, m
        cdef int i,ns
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_i]
        R=sp.empty([m*self.size,ns])
        libGP.draw(self.s, cint(self.size), cint(ns), X_i.ctypes.data_as(ctpd), (cint*len(D))(*D),R.ctypes.data_as(ctpd),cint(m))
        return R
    
    def draw_post(self,X_i,D_i,z):
        cdef int ns
        ns = X_i.shape[0]
        [m,V] = self.infer_full_post(X_i,D_i)
        R = sp.empty([ns,z])
        libGP.drawk(V.ctypes.data_as(ctpd),cint(ns),R.ctypes.data_as(ctpd),cint(z))
        R+=sp.hstack([m.T]*z)
        return R.T
    
    def llk(self):<commit_msg>fixed post_full with single kernel<commit_after># Gaussian process with derivatives in c
#LK_only can only be used for log likelihood or posterior log likelihood
#GPcore takes a single hyperparameter value or a set. With a set the infer_ methods produce individual inferences while infer_post methods return the posterior
#cython: profile=True
#

__author__ = "mark"
__date__ = "$22-Nov-2015 21:27:19$"

#TODO draw doesn't work correclty on a GP created with no data
import numpy as np
cimport numpy as np
import scipy as sp
import ctypes as ct
import os
import sys
from copy import deepcopy as dc
from scipy.stats import norm as norms
from libc.math cimport log10, log, isnan
#print os.path.join(os.path.split(__file__)[0],'../../dist/Release/GNU-Linux/libGPshared.so')
from . import __file__ as fl
from scipy import linalg as spl

libGP = ct.cdll.LoadLibrary(os.path.join(os.path.split(fl)[0],'../cproj/libcproj.so')) #path to c-shared library
#print libGP
libGP.k.restype = ct.c_double

ctpd = ct.POINTER(ct.c_double)
cint = ct.c_int


class GP_LKonly:
    def __init__(self, X_s, Y_s, S_s, D_s, kf):
        cdef int n, D, i
        n ,D = X_s.shape
        self.hyp = kf.hyp
        R = ct.c_double()
        Dc = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_s]
        libGP.newGP_LKonly(cint(D),cint(n),X_s.ctypes.data_as(ctpd),Y_s.ctypes.data_as(ctpd),S_s.ctypes.data_as(ctpd),(cint*len(Dc))(*Dc), cint(kf.Kindex), kf.hyp.ctypes.data_as(ctpd),ct.byref(R))
        self.l = R.value
        
        return
    
    def llk(self):
        #log likelihood of data given hyperparametrs
        return self.l
    
    def plk(self,lm,ls):
        #log likelihood given lognormal prior over hyperparameters
        tmp = 0.

        for i,h in enumerate(self.hyp):
            tmp -= 0.5*((log10(h)-lm[i])**2)/ls[i]**2
        
        return self.l+tmp

class GPcore:
    def __init__(self, X_s, Y_s, S_s, D_s, kf):
        #print [X_s, Y_s, S_s, D_s, kf]
        if type(kf) is kernel:
            self.size = 1
            kf = [kf]
        else:
            self.size = len(kf)
        allhyp = sp.hstack([k.hyp for k in kf])
        self.kf=kf
        [self.n ,self.D] = X_s.shape
        Dx = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_s]
        self.s = libGP.newGP_hypset(cint(self.D),cint(self.n),cint(kf[0].Kindex),X_s.ctypes.data_as(ctpd),Y_s.ctypes.data_as(ctpd),S_s.ctypes.data_as(ctpd),(cint*len(Dx))(*Dx),allhyp.ctypes.data_as(ctpd),cint(self.size))
        self.Y_s=Y_s
        
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_s]
        #print self.get_cho()
        libGP.presolv(self.s,cint(self.size))
        
        return
    
    def __del__(self):
        #print "here"
        libGP.killGP(self.s,cint(self.size))
        #print "also here"
        return
    
    def printc(self):
        print self.size
        libGP.ping(self.s, cint(self.size))
        return
    def get_cho(self):
        C=sp.empty([self.n,self.n*self.size])
        libGP.get_cho(self.s,cint(self.size),C.ctypes.data_as(ctpd))
        return C
    def infer_m(self,X_i,D_i):
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_i]
        R=sp.vstack([sp.empty(ns)]*self.size)
        libGP.infer_m(self.s, cint(self.size), ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
        return R
    
    def infer_m_partial(self,X_i,D_i,ki,hyp):
        
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_i]
        R=sp.vstack([sp.empty(ns)]*1)
        
        #libGP.infer_m_partial(self.s, cint(ki),hyp.ctypes.data_as(ctpd),ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
        libGP.infer_m_partial(self.s,cint(ki),hyp.ctypes.data_as(ctpd),ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
            
        
        return R
    
    def infer_m_post(self,X_i,D_i):
        ns=X_i.shape[0]
        R = self.infer_m(X_i,D_i)
        
        return sp.mean(R,axis=0).reshape([1,ns])
    
    
    def infer_full(self,X_i,D_i):
        cdef int ns,j,i
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**j for j in x])) for x in D_i]
        R=sp.vstack([sp.empty([ns+1,ns])]*self.size)
        libGP.infer_full(self.s, cint(self.size), ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))
        m = sp.vstack([R[i*(ns+1),:] for i in range(self.size)])
        V = sp.vstack([R[(ns+1)*i+1:(ns+1)*(i+1),:] for i in range(self.size)])
        return [m,V]
    
    def infer_full_post(self,X_i,D_i):
        class MJMError(Exception):
            pass
        [m,V] = self.infer_full(X_i,D_i)
        cdef int i,ns
        ns=X_i.shape[0]
        cv = sp.zeros([ns,ns])

        for i in range(self.size):
            cv+=V[ns*i:ns*(i+1),:]
        cv= cv/self.size
        if self.size>1:
            cv+=sp.cov(m,rowvar=0,bias=1)
        #print "_________________"
        #print self.size
        #print sp.cov(m,rowvar=0,bias=1)
        #print V
        #print cv
        return [sp.mean(m,axis=0).reshape([1,ns]),cv]
    
    def infer_diag(self,X_i,D_i):
        cdef int i,j,ns
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**j for j in x])) for x in D_i]
        R=sp.vstack([sp.empty([2,ns])]*self.size)
        libGP.infer_diag(self.s,cint(self.size), ns,X_i.ctypes.data_as(ctpd),(cint*len(D))(*D),R.ctypes.data_as(ctpd))

        m = sp.vstack([R[i*2,:] for i in range(self.size)])
        V = sp.vstack([R[i*2+1,:] for i in range(self.size)])
        return [m,V]
    
    def infer_diag_post(self,X_ii,D_i):
        cdef int ns
        X_i = dc(X_ii)
        ns = len(D_i)
        
        X_i.resize([ns,self.D])
        [m,V] = self.infer_diag(X_i,D_i)
        if sp.amin(V)<=-0.:
            class MJMError(Exception):
                pass
            print "negative/eq variance"
            print [m,V,X_i,D_i]
            print "_______________"
            #self.printc()
            raise(MJMError)
        if sp.amin(sp.var(m,axis=0))<-0.:
            class MJMError(Exception):
                pass
            print "negativevar of mean"
            print X_i.shape
            print [m,V,sp.var(m,axis=0),X_i,D_i]
            print "_______________"
            #self.printc()
            raise(MJMError)
        
        return [sp.mean(m,axis=0).reshape([1,ns]),(sp.mean(V,axis=0)+sp.var(m,axis=0)).reshape([1,ns])]
        
    
    def draw(self,X_i,D_i,m):
        #make m draws at X_i Nd, X, D, R, m
        cdef int i,ns
        ns=X_i.shape[0]
        D = [0 if isnan(x[0]) else int(sum([8**i for i in x])) for x in D_i]
        R=sp.empty([m*self.size,ns])
        libGP.draw(self.s, cint(self.size), cint(ns), X_i.ctypes.data_as(ctpd), (cint*len(D))(*D),R.ctypes.data_as(ctpd),cint(m))
        return R
    
    def draw_post(self,X_i,D_i,z):
        cdef int ns
        ns = X_i.shape[0]
        [m,V] = self.infer_full_post(X_i,D_i)
        R = sp.empty([ns,z])
        libGP.drawk(V.ctypes.data_as(ctpd),cint(ns),R.ctypes.data_as(ctpd),cint(z))
        R+=sp.hstack([m.T]*z)
     #R=sp.random.multivariate_normal(m.flatten(),V,z)
        return R.T
    
    def llk(self):
{'fraction_non_alphanumeric': 0.16108976597974153, 'fraction_numerical': 0.011945511701012924, 'mean_word_length': 2.3425169273873454, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2051337', 'n_tokens_mistral': 6036, 'n_tokens_neox': 5930, 'n_words': 1087}