ordered.src

没有说明

/*
** ordered.src - Ordered Logit and Probit Analysis
**
**
** (C) Copyright 1988-1995  Aptech Systems, Inc.
** All Rights Reserved.
**
** This Software Product is PROPRIETARY SOURCE CODE OF APTECH
** SYSTEMS, INC.    This File Header must accompany all files using
** any portion, in whole or in part, of this Source Code.   In
** addition, the right to create such files is strictly limited by
** Section 2.A. of the GAUSS Applications License Agreement
** accompanying this Software Product.
**
** If you wish to distribute any portion of the proprietary Source
** Code, in whole or in part, you must first obtain written
** permission from Aptech Systems.
**
**-------------------**------------------**-------------------**-----------**
**-------------------**------------------**-------------------**-----------**
**
**
** Purpose:  To estimate the ordered probit or logit model using a GAUSS data
**           set.  By default the ordered probit model is estimated.  The
**           ordered logit model is estimated by setting _QRLOGIT to 1.
**
** Format:   { vnames,b,vc,ndtran,pct,meanx,sdx,fit,df,tol }
**             = ORDERED(dataset,depvar,indvars);
**
** Input:    dataset -- string, name of data file
**
**           depvar -- string, name of dependent variable
**                        - or -
**                     scalar, index of dependent variable.
**
**                     The value of depvar will be truncated before analysis.
**                     Thus, 1.4 is treated as category 1.
**
**           indvars -- Kx1 character vector, names of independent variables.
**                         - or -
**                      Kx1 numeric vector, indices of independent variables.
**
**                      The program adds one variable for the constant term.
**
**           Defaults are provided for the following global input
**           variables.  They can be ignored unless you need control
**           over the other options provided by this procedure.
**
**           WARNING:  If you change the defaults in a command file, the new
**           values will apply in the next program you run using ORDERED unless
**           you change them back.  This can be done by running QUANTSET.
**
**           __altnam -- provides alternative names for the variables.
**
**                      if 0 (default), the original names of the variables
**                      are used.
**
**                      if a ((1+NIVAR)x1) character vector, the first name in
**                      this vector will be used to label the dependent
**                      variable and the remaining NIVAR names will be used to
**                      label the independent variables.
**
**           __miss -- global scalar, default 1.
**
**                     if 0, there are no missing values (fastest).
**
**                     if 1, do listwise deletion, drop an observation
**                     if there are any missing values among the independent
**                     and dependent variables.
**
**           __output -- global scalar, default 1.
**
**                       if 1, sends results to the output device (including
**                       the screen).
**
**                       if 0, no information is sent to output.
**
**           __range  -- 2*1 vector.  The range of record in data
**                       set used for analysis.  The first element is the
**                       starting row index, the second element is the
**                       endding row index.
**
**                       Default is the whole dataset.
**
**           __row -- global scalar, default 0.
**
**                    if 0, the number of rows to read per iteration of the
**                    read loop is calculated by the program.
**
**                    if not 0, the specified number of rows will be read.
**
**           _dtsel -- global scalar, default 0.
**
**                       if 0, all cases are selected for analysis.
**
**                       if Kx3, cases are selected into samples according
**                       to specified conditions.  See DTRAN for details.
**
**           __tol -- global scalar controlling the iterations.  __tol
**                    indicates the maximum difference between estimates
**                    of the coefficients in two adjacent iterations.
**
**          _qrcatnm -- NCATx1 character vector of names of outcome categories
**                       - or -
**                      default scalar 0 in which case names CAT1, CAT2, ...
**                      are used.
**
**          _qrfit -- global scalar, default 0.
**
**                    if 1, print detailed goodness of fit measures, including
**                    table of observed and predicted outcomes.
**
**                    if 0, only print chi-square, -2*log-likehood and percent
**                    correctly predicted.
**
**          _qriter -- global scalar, default 0.
**
**                     if 0, do not print information on iterations.
**
**                     if 1, send detailed information on iterations to
**                     the screen but not to the output device.
**
**                     if 2, send detailed information on iterations to
**                     the output device.
**
**          _qrlogit -- global scalar, default 0;
**
**                      if 1, the ordered logit model is estimated;
**
**                      if 0, the ordered probit model is estimated.
**
**          _qrpred -- global scalar, default 0.
**
**                     if 0, predicted values will not be written to disk.
**
**                     if not 0, predicted probabilities for each outcome
**                     category are written to file ^_qrpred with
**                     NCAT+1 variables.  The first ncat are PRED1,PRED2,
**                     ...,PREDNCAT.  The last variable is the variable
**                     defined by the variable depvar.
**
**          _qrpredn -- string name of dataset for predicted values.  The
**                      default name is "_qrpred".
**
**          _qrstart -- global scalar, default 0.
**
**                      if 0, do not use user supplied start values.
**
**                      if not 0, user should provide a (NCAT-1+NIVAR)
**                       vector of start values.  First, provide start
**                       values for the intercepts, then the slopes.
**
**          _qrstat -- global scalar, default 0.
**
**                     if 0, do not print descriptive statistics.
**
**                     if 1, print descriptive statistics.
**
**          ORDERED uses the method of scoring for estimation, with
**          squeezes.  Squeezes are controlled with these globals:
**
**          _qrsqtol -- global scalar, default .01
**
**                      when the proportional change in the likelihood
**                      function is smaller than _qrsqtol or the change in
**                      the likelihood function is in the wrong direction,
**                      take a squeeze.
**
**          _qrnsqz0 -- global scalar, default 0.
**
**                      if 0, squeezes will not be computed until changes
**                      in the likelihood function from one iteration to
**                      the next become small.
**
**                      if not 0, the program will take up to that number of
**                      squeezes per iteration starting with the first
**                      iteration.
**
**                      Since squeezes take time and are less effective
**                      when estimates are far from the converged values, it
**                      is generally best to leave this as 0.
**
**          _qrsqz -- global scalar, default 0.
**
**                    if 0, don't take squeezes until the change in the
**                    likelihood function is small.
**
**                    if 1, consider taking squeezes from the first iteration.
**
**          _qrnsqz1 -- global scalar, default 10.
**
**                      when squeezes begin, this is the maximum number of
**                      squeezes that will be taken before proceeding to the
**                      next iteration.
**
**          _qrmiter -- maximum number of iterations, default = 1000.
**
**
** Output:   vnames -- a (K+2)x1 character vector containing the names of
**                     the variables in the model.  The order is:
**                     depvar|"CONSTANT"|indvars.
**
**           b -- an NPARM=(NCAT-1)*(K+1) vector of parameter estimates in
**                the order:  intercepts|var1|var2|...varK.  For each
**                variable the parameters are in the order comparing the
**                first category to NCAT, the second to NCAT, ... to
**                NCAT-1 to NCAT.  See below for details.
**
**                If errors are encountered a message will be sent to the
**                error log.  Also, b will contain a scalar error code.  This
**                code appears as missing unless it is translated with
**                the command scalerr(b).  The codes are defined as:
**
**                      1  data file not found
**                      2  found undefined variables
**                     30  system singular
**                     31  too few nonmissing observations.
**                     71  number of categories of dependent variable is less
**                         than 2
**                     72  one of the outcome categories has no cases
**                     73  an independent variable has no variation
**                     74  can't open file for predicted values
**                     75  out of disk space
**                     77  all cases were deleted
**                     78  singular matrix encountered during iterations
**                     79  wrong number of start values specified
**
**
**           vc -- NPARMxNPARM variance covariance matrix for the parameters
**                 in b.
**
**           ndtran -- 2x1 vector of observations.  Element 1 contains
**                     number of cases read from dataset; element 2
**                     contains number of cases left after deletion of
**                     missing cases controlled by __miss, it is the
**                     number of cases used in the analysis.
**
**           pct -- the percent of cases in each of the outcome categories.
**                  Arranged in order lowest to highest.
**
**           meanx -- the means based on nused cases of the independent
**                    variables in the order in indvars.
**
**           sdx -- the standard deviations based on nused cases of the
**                  independent variables in the order in indvars.
**
**           fit -- 4x1 vector of goodness of fit measures.  Element 1 is
**                  the likelihood ratio chi-square assessing the overall
**                  fit of the model; element 2 is -2 times the log
**                  likelihood function evaluated at the estimated values;
**                  element 3 is -2 times the log likelihood function
**                  evaluated with the slopes fixed to zero; element 4 is
**                  the percentage of correct predictions from the model.
**
**           df -- the degrees of freedom associated with lrx2.
**
**           tol -- the tolerance reached.  If convergence was obtained,
**                  tol must be less than __tol.
**
** Remarks:  See the manual for details on the model.
**
** Library:  QUANTAL
**
** See Also: LOGIT, PROBIT, DTRAN
*/

#include gauss.ext
#include quantal.ext

proc(10) = ordered(dataset,depvar,indvars);

    local astrt,abadj,abchng,abml,abnew,abold,absz, absqz,abstrt,bstrt,bm,
        cab,ccdfd,ccdfl,ccdfu,cll,cxb, cdfd,cd,dash,df,dta,start, err,
        errmsg,fin, fl,fmt,fpred,g,i,ii,inf, iter,iters,i1,i2, lbl,llchng,
        llnew,llold,lls,llwrk, m,mask,mask1,mask2,lastobs, maxx,meanx,minx,
        n,ncon,range,nout,nparm,nr,nused,nvar,obspred, omat, p_ml,pct,pcfd,
        pcfl,pcfu,pr,pu,pl,prednm,s, sdx,seml,spc, success,sqz,stplnth,
        sttime,tl,tml,tu, tol,vcml,x,xb,xx,xy,y, readdisk,nm,ndtran,
        oldtrap,xlbl,ylbl,ncat,llrest,llfull,nb, lrx2,lrxp, cdfu,cdfl,ctu,
        ctl,abs1,abs2,absadj,ycat,inm,wx, cx,ydum,oldnfmt,oldcfmt,iswt,
        wlbl,windx,wt,wnused,varindx,nd, count, counter, ky;

    oldnfmt = formatnv("");
    oldcfmt = formatcv("");
    fin = -1;

    sttime = hsec;
    inf = 30;       /* Infinity */
    ndtran = zeros(2,1);
    clear readdisk,nused, wnused;
    dataset = "" $+ dataset;
    open fin = ^dataset;
    if fin == -1;
        errmsg = "ERROR:  Can't open file " $+ dataset $+ ".";
        goto errout(error(1));
    endif;

    { nm, varindx } = indices(dataset, depvar|indvars);
    if scalerr(nm);
        goto errout(nm);
    endif;
    ylbl = nm[1];
    xlbl = trimr(nm,1,0);
    nvar = rows(xlbl);
    if __weight $== 0;
        wlbl = "";
        windx = 0;
        iswt = 0;
    else;           /* INCLUDES __WEIGHT IF DEFINED */
        { wlbl, windx } = indices(dataset, __weight);
        if scalerr(wlbl);
            goto errout(wlbl);
        endif;
        iswt = 1;
        nm = nm|wlbl;
        varindx = varindx|windx;
    endif;
    _qrycat = miss(0,0);

    { nr,start,counter,lastobs } = _rngchk(dataset,__range);
    range = lastobs - start + 1;
    ndtran = zeros(2,1);
    ndtran[1] = range;
    call seekr(fin,start);
    nd = 0;
    count = counter;
    do while count < lastobs;
        readdisk = readdisk+1;
        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;
        nd = rows(dta);
        if nd == 0;
            continue;
        endif;
        ndtran[2] = ndtran[2] + nd;
        y = trunc(dta[.,1]);        /* number of outcome category */
        if iswt;
            wt = dta[.,cols(dta)];
        else;
            wt = ones(rows(dta),1);
        endif;
        wnused = wnused+sumc(wt);
        if scalmiss(_qrycat);
            _qrycat = unique(y,1);
        else;
            _qrycat = union(y,_qrycat,1);
        endif;
    endo;

    nused = ndtran[2];
    if nused == 0;
        errmsg = "ERROR:  No observations left after deleting missing values.";
        goto errout(error(77));
    endif;
    readdisk = readdisk>1;
    ncat = rows(_qrycat);
    ycat = _qrycat[1:ncat]';
    if ncat le 1;
        errmsg = "ERROR:  Too few dependent categories.";
        goto errout(error(71));
    endif;
    if rows(_qrcatnm) .ne ncat;
        _qrtmp = 0 $+ "Group "$+ftocv(_qrycat,1,0);
        inm = 0;
    else;
        _qrtmp = _qrcatnm;
        inm = 1;