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

% Calculates negative loglikelihood for ordinal prbit model

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));
        
    elseif y(i) == J
        z2       = allgammaCp(J)-mu(i);
        lnpdf(i) = log(1-normcdf(z2,0,1));
    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));
    end
end

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

        