function retRes = ginv(inEta,inLink)
%function retRes = ginv(inEta,inLink)
% 
% Calculate the value of g^{-1}(eta) for a generalized linear model with
% link function g(.),
%       g(E[Y}) = mu + X*beta;
%
% Input:
%   inEta: eta=mu+X*beta;
%   inLink: the link function g(.).
%
% Output:
%   retRes: g^{-1}(eta).
%
% Dabao Zhang    May 8,2013;
%
% Revised by:
%       Min Ren, March 17, 2015 -- Implement Multinomial proportinal odds &
%                                  continuation ratio models
%
% Copyright (c) 2013 by Dabao Zhang.

error(nargchk(2, 2, nargin));

switch inLink
    case 'logit'
        % -36<eta<36 as we want eps<1/(1+exp(-eta))<1-eps
        inEta(inEta>36) = 36;
        inEta(inEta<-36) = -36;
        retRes = 1./(1+exp(-inEta));

    case 'probit'
        retRes = normcdf(inEta,0,1);

    case 'comploglog'
        % -350<eta<350 as we want realmin<exp(eta*2))<realmax
        inEta(inEta>350) = 350;
        inEta(inEta<-350) = -350;
        retRes = 1-exp(-exp(inEta));

    case 'identity'
        retRes = inEta;

    case 'log'
        % -350<eta<350 as we want realmin<exp(eta*2))<realmax
        inEta(inEta>350) = 350;
        inEta(inEta<-350) = -350;
        retRes = exp(inEta);

    case 'loglog'
        % realmin<exp(eta*2))<exp(2*exp(eta))<realmax
        inEta(inEta>5.8) = 5.8;
        inEta(inEta<-350) = -350;
        retRes = exp(exp(inEta));
    
    case 'bclogit'
        maxEta = max(inEta,[],2);
        maxEta(maxEta<0) = 0;
        retRes = exp(inEta-maxEta);
        denom = exp(-maxEta)+sum(retRes,2);
        retRes = retRes/denom;
        
    case 'pologit'
        % -36<eta<36 as we want eps<1/(1+exp(-eta))<1-eps
        %inEta(inEta>36) = 36;
        %inEta(inEta<-36) = -36;
        p = 1./(1+exp(-inEta));
        % p_1,p_1,...,p_k (p_0 is excluded)
        retRes = [p(:,2:end)-p(:,1:end-1) 1-p(:,end)];
    
    case 'crlogit'
        p = 1./(1+exp(inEta));
        picumsum = cumprod(p,2);
        retRes = [picumsum(:,1:(end-1))-picumsum(:,2:end) picumsum(:,end)];

    otherwise
        if( isnumeric(inLink) )
            % eta^(1/k) with k=inLink
            % realmin < (eta^(2/k) < realmax
            minb = (realmin)^(inLink/2);
            maxb = (realmax)^(inLink/2);
            if( inLink<0 )
                maxb = max(maxb,realmin);
                minb = min(minb,realmax)
                inEta(inEta<maxb) = maxb;
                inEta(inEta>minb) = minb;
            elseif( inLink==0 )
                error('    gPOCRE: The link is not supported...');
            elseif( inLink==1 )
                error('    gPOCRE: Please use identity link...');
            else
                minb = max(minb,realmin);
                maxb = min(maxb,realmax);
                inEta(inEta<minb) = minb;
                inEta(inEta>maxb) = maxb;
            end
            retRes = inEta.^(1/inLink);
        else
            error('  gPOCRE: The link is not supported...');
        end
end
