{"action":"created","issue":{"url":"https://api.github.com/repos/GPflow/GPflow/issues/37","repository_url":"https://api.github.com/repos/GPflow/GPflow","labels_url":"https://api.github.com/repos/GPflow/GPflow/issues/37/labels{/name}","comments_url":"https://api.github.com/repos/GPflow/GPflow/issues/37/comments","events_url":"https://api.github.com/repos/GPflow/GPflow/issues/37/events","html_url":"https://github.com/GPflow/GPflow/pull/37","id":148015137,"number":37,"title":"FITC method for GPflow","user":{"login":"alexggmatthews","id":6596998,"avatar_url":"https://avatars.githubusercontent.com/u/6596998?v=3","gravatar_id":"","url":"https://api.github.com/users/alexggmatthews","html_url":"https://github.com/alexggmatthews","followers_url":"https://api.github.com/users/alexggmatthews/followers","following_url":"https://api.github.com/users/alexggmatthews/following{/other_user}","gists_url":"https://api.github.com/users/alexggmatthews/gists{/gist_id}","starred_url":"https://api.github.com/users/alexggmatthews/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/alexggmatthews/subscriptions","organizations_url":"https://api.github.com/users/alexggmatthews/orgs","repos_url":"https://api.github.com/users/alexggmatthews/repos","events_url":"https://api.github.com/users/alexggmatthews/events{/privacy}","received_events_url":"https://api.github.com/users/alexggmatthews/received_events","type":"User","site_admin":false},"labels":[{"url":"https://api.github.com/repos/GPflow/GPflow/labels/help%20wanted","name":"help wanted","color":"159818"}],"state":"open","locked":false,"assignee":{"login":"jameshensman","id":154383,"avatar_url":"https://avatars.githubusercontent.com/u/154383?v=3","gravatar_id":"","url":"https://api.github.com/users/jameshensman","html_url":"https://github.com/jameshensman","followers_url":"https://api.github.com/users/jameshensman/followers","following_url":"https://api.github.com/users/jameshensman/following{/other_user}","gists_url":"https://api.github.com/users/jameshensman/gists{/gist_id}","starred_url":"https://api.github.com/users/jameshensman/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/jameshensman/subscriptions","organizations_url":"https://api.github.com/users/jameshensman/orgs","repos_url":"https://api.github.com/users/jameshensman/repos","events_url":"https://api.github.com/users/jameshensman/events{/privacy}","received_events_url":"https://api.github.com/users/jameshensman/received_events","type":"User","site_admin":false},"milestone":null,"comments":14,"created_at":"2016-04-13T10:21:43Z","updated_at":"2016-04-13T15:52:05Z","closed_at":null,"pull_request":{"url":"https://api.github.com/repos/GPflow/GPflow/pulls/37","html_url":"https://github.com/GPflow/GPflow/pull/37","diff_url":"https://github.com/GPflow/GPflow/pull/37.diff","patch_url":"https://github.com/GPflow/GPflow/pull/37.patch"},"body":"This is not ready to merge yet. Rather it is a \"help wanted\" PR for @jameshensman to look at.\r\n\r\nI working from Carl and Hannes' GPML code:\r\n\r\n```\r\nfunction [post nlZ dnlZ] = infFITC(hyp, mean, cov, lik, x, y)\r\n\r\n% FITC approximation to the posterior Gaussian process. The function is\r\n% equivalent to infExact with the covariance function:\r\n%   Kt = Q + G; G = diag(g); g = diag(K-Q);  Q = Ku'*inv(Quu)*Ku;\r\n% where Ku and Kuu are covariances w.r.t. to inducing inputs xu, snu2 = sn2/1e6\r\n% is the noise of the inducing inputs and Quu = Kuu + snu2*eye(nu).\r\n% We fixed the standard deviation of the inducing inputs snu to be a one per mil\r\n% of the measurement noise's standard deviation sn.\r\n% The implementation exploits the Woodbury matrix identity\r\n%   inv(Kt) = inv(G) - inv(G)*V'*inv(eye(nu)+V*inv(G)*V')*V*inv(G)\r\n% in order to be applicable to large datasets. The computational complexity\r\n% is O(n nu^2) where n is the number of data points x and nu the number of\r\n% inducing inputs in xu.\r\n%\r\n% The function takes a specified covariance function (see covFunction.m) and\r\n% likelihood function (see likFunction.m), and is designed to be used with\r\n% gp.m and in conjunction with covFITC and likGauss.\r\n%\r\n% The inducing points can be specified through 1) the 2nd covFITC parameter or\r\n% by 2) providing a hyp.xu hyperparameters. Note that 2) has priority over 1).\r\n% In case 2) is provided and derivatives dnlZ are requested, there will also be\r\n% a dnlZ.xu field allowing to optimise w.r.t. to the inducing points xu. However\r\n% the derivatives dnlZ.xu can only be computed for one of the following eight\r\n% covariance functions: cov{Matern|PP|RQ|SE}{iso|ard}.\r\n%\r\n% Copyright (c) by Ed Snelson, Carl Edward Rasmussen \r\n%                                               and Hannes Nickisch, 2013-10-28.\r\n%\r\n% See also INFMETHODS.M, COVFITC.M.\r\n\r\nif iscell(lik), likstr = lik{1}; else likstr = lik; end\r\nif ~ischar(likstr), likstr = func2str(likstr); end\r\nif ~strcmp(likstr,'likGauss')               % NOTE: no explicit call to likGauss\r\n  error('Inference with inFITC only possible with Gaussian likelihood.');\r\nend\r\ncov1 = cov{1}; if isa(cov1, 'function_handle'), cov1 = func2str(cov1); end\r\nif ~strcmp(cov1,'covFITC'); error('Only covFITC supported.'), end    % check cov\r\nif isfield(hyp,'xu'), cov{3} = hyp.xu; end  % hyp.xu is provided, replace cov{3}\r\n\r\n[diagK,Kuu,Ku] = feval(cov{:}, hyp.cov, x);         % evaluate covariance matrix\r\nm = feval(mean{:}, hyp.mean, x);                          % evaluate mean vector\r\n[n, D] = size(x); nu = size(Kuu,1);\r\ncov{3} = reshape(cov{3},[],D);\r\n\r\nsn2  = exp(2*hyp.lik);                              % noise variance of likGauss\r\nsnu2 = 1e-6*sn2;                              % hard coded inducing inputs noise\r\nLuu  = chol(Kuu+snu2*eye(nu));                         % Kuu + snu2*I = Luu'*Luu\r\nV  = Luu'\\Ku;                                     % V = inv(Luu')*Ku => V'*V = Q\r\ng_sn2 = diagK + sn2 - sum(V.*V,1)';          % g + sn2 = diag(K) + sn2 - diag(Q)\r\nLu = chol(eye(nu) + (V./repmat(g_sn2',nu,1))*V');  % Lu'*Lu=I+V*diag(1/g_sn2)*V'\r\nr  = (y-m)./sqrt(g_sn2);\r\nbe = Lu'\\(V*(r./sqrt(g_sn2)));\r\niKuu = solve_chol(Luu,eye(nu));                       % inv(Kuu + snu2*I) = iKuu\r\npost.alpha = Luu\\(Lu\\be);                      % return the posterior parameters\r\npost.L  = solve_chol(Lu*Luu,eye(nu)) - iKuu;                    % Sigma-inv(Kuu)\r\npost.sW = ones(n,1)/sqrt(sn2);           % unused for FITC prediction  with gp.m\r\n\r\nif nargout>1                                % do we want the marginal likelihood\r\n  nlZ = sum(log(diag(Lu))) + (sum(log(g_sn2)) + n*log(2*pi) + r'*r - be'*be)/2;\r\n  if nargout>2                                         % do we want derivatives?\r\n    dnlZ = hyp;                                 % allocate space for derivatives\r\n    al = r./sqrt(g_sn2) - (V'*(Lu\\be))./g_sn2;      % al = (Kt+sn2*eye(n))\\(y-m)\r\n    B = iKuu*Ku; w = B*al;\r\n    W = Lu'\\(V./repmat(g_sn2',nu,1));\r\n    for i = 1:numel(hyp.cov)\r\n      [ddiagKi,dKuui,dKui] = feval(cov{:}, hyp.cov, x, [], i);  % eval cov deriv\r\n      R = 2*dKui-dKuui*B; v = ddiagKi - sum(R.*B,1)';   % diag part of cov deriv\r\n      dnlZ.cov(i) = (ddiagKi'*(1./g_sn2)+w'*(dKuui*w-2*(dKui*al))-al'*(v.*al)...\r\n                                 - sum(W.*W,1)*v - sum(sum((R*W').*(B*W'))) )/2;\r\n    end\r\n    diag_dK = 1./g_sn2 - sum(W.*W,1)' - al.*al;                  % diag(dnlZ/dK)\r\n    dnlZ.lik = sn2*sum(diag_dK);\r\n    % since snu2 is a fixed fraction of sn2, there is a covariance-like term in\r\n    BWt = B*W'; % the derivative as well\r\n    dKuui = 2*snu2; R = -dKuui*B; v = -sum(R.*B,1)';    % diag part of cov deriv\r\n    dnlZ.lik = dnlZ.lik + (w'*dKuui*w -al'*(v.*al)...\r\n                                    - sum(W.*W,1)*v - sum(sum((R*W').*BWt)) )/2; \r\n    for i = 1:numel(hyp.mean)\r\n      dnlZ.mean(i) = -feval(mean{:}, hyp.mean, x, i)'*al;\r\n    end\r\n    if isfield(hyp,'xu')                 % derivatives w.r.t. inducing points xu\r\n      xu = cov{3};\r\n      cov = cov{2};           % get the non FITC part of the covariance function\r\n      Kpu  = cov_deriv_sq_dist(cov,hyp.cov,xu,x);           % d K(xu,x ) / d D^2\r\n      Kpuu = cov_deriv_sq_dist(cov,hyp.cov,xu);             % d K(xu,xu) / d D^2\r\n      if iscell(cov), covstr = cov{1}; else covstr = cov; end\r\n      if ~ischar(covstr), covstr = func2str(covstr); end\r\n      if numel(strfind(covstr,'iso'))>0            % characteristic length scale\r\n        e = 2*exp(-2*hyp.cov(1));\r\n      else\r\n        e = 2*exp(-2*hyp.cov(1:D));\r\n      end\r\n      v = diag_dK-1./g_sn2;        % BdK = B * ( dnlZ/dK - diag(diag(dnlZ/dK)) )\r\n      BdK = B.*repmat(v',nu,1) + BWt*W + (B*al)*al';\r\n      A = Kpu.*BdK; C = Kpuu.*(BdK*B'); C = diag(sum(C,2)-sum(A,2)) - C;\r\n      dnlZ.xu = A*x*diag(e) + C*xu*diag(e);  % bring in data and inducing points\r\n    end\r\n  end\r\nend\r\n\r\n```\r\n(Sorry to paste apparently .m files aren't allowed as attachments). \r\n\r\nCurrently I get numerical inversion problems when I try to run the code. \r\n\r\nThe tests currently fail. \r\n\r\nAny ideas @jameshensman ?\r\n"},"comment":{"url":"https://api.github.com/repos/GPflow/GPflow/issues/comments/209521009","html_url":"https://github.com/GPflow/GPflow/pull/37#issuecomment-209521009","issue_url":"https://api.github.com/repos/GPflow/GPflow/issues/37","id":209521009,"user":{"login":"alexggmatthews","id":6596998,"avatar_url":"https://avatars.githubusercontent.com/u/6596998?v=3","gravatar_id":"","url":"https://api.github.com/users/alexggmatthews","html_url":"https://github.com/alexggmatthews","followers_url":"https://api.github.com/users/alexggmatthews/followers","following_url":"https://api.github.com/users/alexggmatthews/following{/other_user}","gists_url":"https://api.github.com/users/alexggmatthews/gists{/gist_id}","starred_url":"https://api.github.com/users/alexggmatthews/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/alexggmatthews/subscriptions","organizations_url":"https://api.github.com/users/alexggmatthews/orgs","repos_url":"https://api.github.com/users/alexggmatthews/repos","events_url":"https://api.github.com/users/alexggmatthews/events{/privacy}","received_events_url":"https://api.github.com/users/alexggmatthews/received_events","type":"User","site_admin":false},"created_at":"2016-04-13T15:52:05Z","updated_at":"2016-04-13T15:52:05Z","body":"I just talked to Rich Turner about it. He says do what you would do for variational. i.e integrate out q(u) using the full conditional.\r\n\r\nIncidentally isn't GitHub useful?\r\n\r\n@markvdw @jameshensman  I think I agree with Mark's second point but I'm not sure I think it is that important. I'd rather have q treat train and test points differently and give a better prediction. I'm not sure I agree with Mark's first point. I think prediction goes as O(nHoldOut^2) rather than O(nHoldOut^3)  (neglecting M dependence) . It is still interesting to look at correlation among a small set of hold out points and in any case we do supply a version that just gives the diagonals of the covariance and scales as O(nHoldOut). "}}