
function [negsumLogL] = OPnegSumloglike(delta,y,x,beta)

% Calculates negative loglikelihood for ordinal prbit model
% This function was earlier named as OPnegloglikelihood, when working in
% Rahman(2016).
% Output
%--------------------------------------------------------------------------
% NegSumLnPdf : Negative sum of the log-likelihood
% Input
% delta     : cut-points, a column vector
% y         : observed dicrete values, a column vector
% x         : covariate matrix including a column of ones
% beta      : Gibbs draw of beta, a column vector
%--------------------------------------------------------------------------

n           = size(x,1);
J           = size(unique(y),1);
lnpdf       = zeros(n,1);

if J < 4
    error('Number of outcomes should be ge 4')
end

% Transforming the delta to original cut-points
% delta(j) = ln (gamma(j) - gamma(j-1))

expdelta        = exp(delta);
gammaCp         = zeros(size(expdelta,1)+1,1);

for j = 2:J-1
    gammaCp(j)  = sum(expdelta(1:j-1,1));
end
allgammaCp      = [-Inf; gammaCp; Inf];

mu              = x*beta;

for i = 1:n
    if y(i) == 1
        z1       = -mu(i);
        %lnpdf(i) = log(normcdf(z1,0,1));
        lnpdf(i) = log(fast_normcdf(z1));
        
    elseif y(i) == J
        z2       = allgammaCp(J)-mu(i);
        % lnpdf(i) = log(1-normcdf(z2,0,1));
        lnpdf(i) = log(1-fast_normcdf(z2));

    else
        z3       = allgammaCp(y(i)+1)  - mu(i);
        z4       = allgammaCp(y(i))    - mu(i);
        %lnpdf(i) = log(normcdf(z3,0,1) - normcdf(z4,0,1));
        lnpdf(i) = log(fast_normcdf(z3) - fast_normcdf(z4));
    end
end

negLogL    = - lnpdf;
negsumLogL = sum(negLogL);

end