ordered.src

没有说明

        errorlog "WARNING:  Maximum iterations exceeded";
    elseif ky;
        errorlog "WARNING:  convergence forced";
    endif;

    if _qrpred;

    /* Open file for predicted probabilities */

        if inm;
            prednm = (0 $+ "P_" $+ _qrtmp)|ylbl;
        else;
            prednm = (0 $+ "P_" $+ ftocv(_qrycat,1,0))|ylbl;
        endif;

        create fpred = ^_qrpredn with ^prednm,ncat+1,8;
        if fpred == -1;
            errmsg = "Can't open file for residuals: " $+ _qrpredn $+ ".";
            goto errout(error(74));
        endif;
    endif;

    if readdisk;
        call seekr(fin,start);
        count = counter;
        do while count < lastobs;
            dta = readr(fin,nr);
            count = count + rows(dta);
            if count > lastobs;
                dta = trimr(dta,0,count-lastobs);
            endif;
            if __miss == 1;
                dta = packr(dta[.,varindx]);
            else;
                dta = dta[.,varindx];
            endif;
            if scalmiss(dta);
                continue;
            endif;
            ydum = trunc(dta[.,1]);
            y = maxindc(((ones(rows(dta),1)*ycat .== ydum))');
            if iswt;
                wt = dta[.,cols(dta)];
                dta = trimr(dta',0,1)';
            else;
                wt = ones(rows(dta),1);
            endif;
            x = trimr(dta',1,0)';
            gosub pred;
            if _qrpred;     /* write observed outcome and probabilities  */
                if writer(fpred,pr~ydum) /= rows(ydum);
                    errmsg = "ERROR:  Disk full.  Predicted file not comple"\
                        "ted.";
                    goto errout(error(75));
                endif;
            endif;
        endo;
    else;
        gosub pred;
        if _qrpred;         /* write observed outcome and probabilities  */
            if writer(fpred,pr~ydum) /= rows(ydum);
                errmsg = "ERROR:  Disk full.  Predicted file not completed.";
                goto errout(error(75));
            endif;
        endif;
    endif;
    obspred = reshape(obspred',ncat,ncat)./wnused;
    /* percent correct predictions - Maddall, p.70 */
    success = 100*sumc(diag(obspred));
    obspred = wnused*(obspred~sumc(obspred'));
    if __output;
        print ("  Percent Correctly Predicted:            " $+
            ftos(success,"*.*lf",10,4));
    endif;

    if __output and _qrfit;
        print;
        print "TABLE OF OBSERVED AND PREDICTED OUTCOMES:";
        print;
        dash = "-----------------------------------------------------------"\
            "------------------";
        spc = "                                                            "\
            "                 ";
        spc = strsect(spc,1,(4.5*ncat-4));
        print "           | " spc "Predicted";
        dash = strsect(dash,1,22+9*ncat);
        call formatcv("*.*lf "~8~8);
        call formatnv("*.*lf "~8~0);
        print "  Observed |";;
        call printfmt(_qrtmp'~"Total",0);
        print;
        print ("  " $+ dash);
        ii = 1;
        call formatcv("*.*lf"~10~8);
        do until ii>ncat;
            call printfmt(_qrtmp[ii],0);
            print " |";;
            call printfmt(obspred[ii,.],1);
            print;
            ii = ii+1;
        endo;
        print ("  " $+ dash);
        print "     Total |";;
        call printfmt(sumc(obspred)',1);
        print;
    endif;

    if _qrpred;
        nout = rowsf(fpred);
        fpred = close(fpred);
        if fpred == 0 and __output;
            print;
            print "PREDICTED VALUES SUCCESSFULLY WRITTEN TO DISK:";
            print;
            print "  The file " _qrpredn " was written with ";;
            print ftos(nout,"*.*lf",1,0) " cases.";
            print;
            print "  The following variables are in the file:";
            print;
            print ("    Prob Y=i for i=1 to " $+ ftos(ncat,"*.*lf",1,0) $+
                ":   ");;
            call printfmt(prednm[1:ncat,1]',0);
            print;
            print;
            print ("    Dependent variable (Y):  " $+ ylbl);
        endif;
    endif;

    print;

    goto done;

/* Subroutine to compute gradient matrix. */

grad:

    if n /= rows(mask);
        mask = reshape(seqa(1,1,ncon+2)',n,ncon+2);
    endif;
    mask1 = mask.==y;
    mask2 = mask.==(y+1);
    xb = x*abold[ncat:ncon+nvar,1];
    tu = submat(-inf|abold[1:ncon,1]|inf,y+1,0)-xb;
    tl = submat(-inf|abold[1:ncon,1]|inf,y,0)-xb;
    clear xb;

    if _qrlogit;
        cdfu = 1./(1+exp(-tu));
        cdfl = 1./(1+exp(-tl));
        pcfu = cdfu.*(1-cdfu);
        pcfl = cdfl.*(1-cdfl);
    else;
        pcfu = pdfn(tu);
        pcfl = pdfn(tl);
        cdfu = cdfn(tu);
        cdfl = cdfn(tl);
    endif;

    clear tu,tl;
    cdfd = cdfu-cdfl;
    pcfd = pcfu-pcfl;

    g = (mask1.*(-pcfl./cdfd))+(mask2.*(pcfu./cdfd));
    g = g[.,2:ncon+1]~-(pcfd./cdfd).*x;
    clear pcfd;

    return;

/* Compute Squeezes:  compares minus the log-likelihood at abold+s*abadj,
   and at abold+(0.5*stplnth*abadj), with abadj = the proposed
   adjustment to the parameters, and with step length stplnth initially 1.
   stplnth is halved until minus the log-likelihood stops declining. */

squeeze:

    pop absadj;
    pop absz;

    abs1 = absz + absadj;
    stplnth = 1/2;
    if minc(abs(absadj)) > _qrsqtol;
        stplnth = meanc(rndu(1,1).*(abs(absadj)));
    endif;

    iters = 1;
    gosub compl(abs1);
    pop llwrk;

    do until iters > _qrnsqz0;
        abs2 = absz+stplnth*absadj;
        gosub compl(abs2);
        pop lls;
        if llwrk >= lls;
            abs1 = abs2;
            stplnth = stplnth/2;
            if minc(abs(absadj)) < _qrsqtol;
                stplnth = stplnth/4;
            endif;
            abs2 = absz + stplnth*absadj;
            llwrk = lls;
            iters = iters+1;
        else;
            sqz = iters-1;
            return(abs1);
        endif;
    endo;
    sqz = 10;
    return(abs2);

/* Subroutine to compute value of log likelihood for squeezes. */

compl:

    pop cab;
    cll = 0;

    if not readdisk;
        cxb = x*cab[ncat:ncon+nvar,1];
        if _qrlogit;
            ctu = submat(-inf|cab[1:ncon,1]|inf,y+1,0)-cxb;
            ctl = submat(-inf|cab[1:ncon,1]|inf,y,0)-cxb;
            ccdfu = 1./(1+exp(-ctu));
            ccdfl = 1./(1+exp(-ctl));
            clear ctu,ctl;
        else;
            ccdfu = cdfn(submat(-inf|cab[1:ncon,1]|inf,y+1,0)-cxb);
            ccdfl = cdfn(submat(-inf|cab[1:ncon,1]|inf,y,0)-cxb);
        endif;
        ccdfd = ccdfu-ccdfl;
        if ccdfd > 0;
            cll = -sumc(wt.*ln(ccdfd));
        else;
            cll = -sumc(wt.*ln(ccdfd-((ccdfd.<=0).*(ccdfd-.0001))));
    /* fix up zero's */
        endif;
    else;
        call seekr(fin,start);
        count = counter;
        do while count < lastobs;
            dta = readr(fin,nr);
            count = count + rows(dta);
            if count > lastobs;
                dta = trimr(dta,0,count-lastobs);
            endif;
            if __miss == 1;
                dta = packr(dta[.,varindx]);
            else;
                dta = dta[.,varindx];
            endif;
            if scalmiss(dta);
                continue;           /* no observations returned */
            endif;
            y = maxindc(((ones(rows(dta),1)*ycat .== trunc(dta[.,1])))');
            if iswt;
                wt = dta[.,cols(dta)];
                dta = trimr(dta',0,1)';
            else;
                wt = ones(rows(dta),1);
            endif;
            x = trimr(dta',1,0)';
            cxb = x*cab[ncat:ncon+nvar,1];
            if _qrlogit;
                ctu = submat(-inf|cab[1:ncon,1]|inf,y+1,0)-cxb;
                ctl = submat(-inf|cab[1:ncon,1]|inf,y,0)-cxb;
                ccdfu = 1./(1+exp(-ctu));
                ccdfl = 1./(1+exp(-ctl));
                clear ctu,ctl;
            else;
                ccdfu = cdfn(submat(-inf|cab[1:ncon,1]|inf,y+1,0)-cxb);
                ccdfl = cdfn(submat(-inf|cab[1:ncon,1]|inf,y,0)-cxb);
            endif;
            ccdfd = ccdfu-ccdfl;
            if ccdfd > 0;
                cll = -sumc(wt.*ln(ccdfd))+cll;
            else;
                cll = -sumc(wt.*ln(ccdfd-((ccdfd.<=0).*(ccdfd-.0001))))+cll;
    /* fix up zero's */
            endif;
        endo;
        i = seekr(fin,1);
    endif;
    return(cll);

pred:

    /* Compute predicted probabilities */

    fl = x*abml[ncon+1:nparm];
    if _qrlogit;
        pr = exp(abml[1]-fl)./(1+exp(abml[1]-fl));
    else;
        pr = cdfn(abml[1]-fl);
    endif;
    i = 2;
    tl = pr;
    do until i .ge ncat;
        if _qrlogit;
            pr = pr~((exp(abml[i]-fl)./(1+exp(abml[i]-fl)))-tl);
        else;
            pr = pr~(cdfn(abml[i]-fl)-tl);
        endif;
        tl = tl+pr[.,i];
        i = i+1;
    endo;

    pr = pr~(ones(rows(x),1) - sumc(pr'));
    /* compute obs/predicted table */
    ii = ((vec((seqa(1,1,ncat)*ones(1,ncat))')*10)+reshape(seqa(1,1,ncat),
        ncat*ncat,1));
    obspred = obspred+countwts(10*y + maxindc(pr'),ii,wt);
    return;

done:

    lrx2 = lrx2|llfull|llrest|success;
    call formatcv(oldcfmt);
    call formatnv(oldnfmt);
    if fin > 0;
        fin = close(fin);
    endif;
    ndpclex;
    retp(ylbl|"CONSTANT"|xlbl,abold,vcml,ndtran,pct,meanx,sdx,lrx2,df,tol);

errout:

    pop err;
    cls;
    if fin > 0;
        fin = close(fin);
    endif;
    ndpclex;
    if not trapchk(1);
        errorlog errmsg;
        print;
        end;
    endif;
    retp(0,err,0,0,0,0,0,0,0,99999);

endp;       /* End of ORDERED */

proc(0) = orditer(iter,sqz,st,en,lbl);
    local omat,fmt,mask,st0,i;
    if sqz ne 0;
        if sqz eq 1;
            print ("ITERATION " $+ ftos(iter,"*.*lf",1,0) $+ " WITH " $+
                ftos(sqz,"*.*lf",1,0) $+ " SQUEEZE:");
        else;
            print ("ITERATION " $+ ftos(iter,"*.*lf",1,0) $+ " WITH " $+
                ftos(sqz,"*.*lf",1,0) $+ " SQUEEZES:");
        endif;
    else;
        print "ITERATION " ftos(iter,"*.*lf",1,0) ":";
    endif;
    print;
    print "  Parameter    Start Value        End Value    % Change";
    print "  -----------------------------------------------------";
    mask = 0~1~1~1;
    let fmt[4,3] = "*.*lf" 10 8 "*.*lf " 16 8 "*.*lf " 16 8 "*.*lf " 11 5;
    i = 0;
    st0 = st;
    do until i == rows(st);
        i = i+1;
        if st[i]/=0;
            st0[i] = 100*(en[i]-st[i])/st[i];
        else;
            st0[i] = -999;
        endif;
    endo;
    omat = lbl~st~en~st0;
    call printfm(omat,mask,fmt);
endp;       /* orditer */