Chapter 2: Regression

RSUBM

Retrieves a symmetric submatrix from a symmetric matrix.

Required Arguments

A — NA by NA symmetric matrix.   (Input)
Only the upper triangle of A is referenced.

SWEPT — Vector of length NA.   (Input)
Element A(I, J) is included in submatrix ASUB if and only if SWEPT(I) > 0.0 and SWEPT(J) > 0.0.

NASUB — Order of submatrix ASUB.   (Output)
NASUB equals the number of elements in SWEPT that are greater than zero.

ASUB — NASUB by NASUB symmetric matrix containing a submatrix of A.   (Output)
If A is not needed, ASUB and A can share the same storage locations.

Optional Arguments

NA — Order of matrix A.   (Input)
Default: NA = size (A,2).

LDA — Leading dimension of A exactly as specified in the dimension statement of the calling program.   (Input)
Default: LDA = size (A,1).

LDASUB — Leading dimension of ASUB exactly as specified in the dimension statement of the calling program.   (Input)
Default: LDASUB = size (ASUB,1).

FORTRAN 90 Interface

Generic:                              CALL RSUBM (A, SWEPT, NASUB, ASUB [,…])

Specific:                             The specific interface names are S_RSUBM and D_RSUBM.

FORTRAN 77 Interface

Single:                                CALL RSUBM (NA, A, LDA, SWEPT, NASUB, ASUB, LDASUB)

Double:                              The double precision name is DRSUBM.

Description

Routine RSUBM retrieves a symmetric submatrix from a symmetric matrix A. If elements i and j of the input vector SWEPT are greater than zero, then the ij-th element of A is output in the submatrix ASUB. Otherwise, the ij-th element of A will not be included in ASUB. (Here, i = 1, 2, …, NA, and j = 1, 2, …, NA, where NA is the order of A.)

Routine RSUBM can be useful in conjunction with two routines, GSWEP and RSTEP. The routine RSUBM can be used after routine GSWEP in order to retrieve the submatrix of A that corresponds to the rows/columns that have been successfully swept. In this case, the SWEPT vector output from GSWEP can be used as the input for the argument SWEPT in RSUBM. Also, RSUBM can be used after routine RSTEP in order to retrieve the submatrix of COVS that corresponds to the independent variables in the final model. In this case, the HIST vector output from RSTEP can be used as the input for the argument SWEPT in RSUBM.

Comments

1.         Workspace may be explicitly provided, if desired, by use of R2UBM/DR2UBM. The reference is:

CALL R2UBM (NA, A, LDA, SWEPT, NASUB, ASUB, LDASUB, IWK)

The additional argument is:

IWK — Vector of length NASUB.

2.         Routine RSUBM can be used after invoking routines GSWEP and RSTEP in order to retrieve the submatrix for the variables in the model.

Example 1

The 2 2 symmetric submatrix ASUB is retrieved from rows and columns 1 and 4 of the 4 4 symmetric matrix A.

 

      USE RSUBM_INT

      USE WRRRN_INT

 

      IMPLICIT   NONE

      INTEGER    LDA, LDASUB, NA

      PARAMETER  (LDASUB=2, NA=4, LDA=NA)

!

      INTEGER    NASUB

      REAL       A(LDA,NA), ASUB(LDASUB,LDASUB), SWEPT(NA)

!

      DATA SWEPT/1.0, -1.0, -1.0, 1.0/

      DATA A/10.0, 20.0, 40.0, 70.0, 20.0, 30.0, 50.0, 80.0, 40.0,&

          50.0, 60.0, 90.0, 70.0, 80.0, 90.0, 100.0/

!

      CALL RSUBM (A, SWEPT, NASUB, ASUB)

      CALL WRRRN ('ASUB', ASUB)

      END

Output

 

      ASUB
        1       2
1    10.0    70.0
2    70.0   100.0

Additional Example

Example 2

This example invokes RSUBM after routine RSTEP in order to retrieve the submatrix of COVS that corresponds to the independent variables in the final stepwise model. With this submatrix, routine BLINF (IMSL MATH/LIBRARY) is used to compute the estimated standard deviation for the intercept in the final model.

A data set from Draper and Smith (1981, pages 629630) is used. The means and the corrected sum of squares and crossproducts matrix for this data are given in the DATA statements. They can be computed using routine CORVC in Chapter 3, “Correlation”. The first four entries in XMEAN and the first four columns of COV correspond to the independent variables, the last entry in XMEAN and the last column of COV correspond to the dependent variable.

After RSTEP is invoked to obtain a model, the intercept is computed using the formula

where k is the number of independent variables in the final model. The estimated standard deviation of the intercept is computed using the formula

where s2 is the error mean square from the fit (stored in AOV(8)), n is the number of observations,  is the subvector of means for the independent variables in the final model (in this case the first mean and the fourth mean), and A is the submatrix (in this case with rows and columns 1 and 4) of the matrix COVS that is output by RSTEP.

 

      USE IMSL_LIBRARIES

 

      IMPLICIT   NONE

      INTEGER    LDCOEF, LDCOV, LDCOVS, NVAR

      PARAMETER  (NVAR=5, LDCOEF=NVAR, LDCOV=NVAR, LDCOVS=NVAR)

!

      INTEGER    I, IEND, INVOKE, IPRINT, ISTEP, J, LEVEL(NVAR), &

                 NFORCE, NIND, NOBS, NOUT, NSTEP

      REAL       AOV(13), B0, COEF(LDCOEF,5), &

                 COV(LDCOV,NVAR), COVS(LDCOVS,NVAR), HIST(NVAR), PIN, &

                 POUT, SCALE(NVAR), SEB0, SQRT, TOL, XMEAN(NVAR)

      INTRINSIC  SQRT

!

      DATA COV/415.231, 251.077, -372.615, -290.000, 775.962, 251.077, &

           2905.69, -166.538, -3041.00, 2292.95, -372.615, -166.538, &

           492.308, 38.0000, -618.231, -290.000, -3041.00, 38.0000, &

           3362.00, -2481.70, 775.962, 2292.95, -618.231, -2481.70, &

           2715.76/

      DATA XMEAN/7.46154, 48.1538, 11.7692, 30.0000, 95.4231/

      DATA LEVEL/4*1, -1/

!

      J = 0

      ISTEP  = 1

      NOBS   = 13

      IPRINT = 1

      CALL RSTEP (COV, NOBS, AOV, COEF, COVS, NVAR=NVAR, ISTEP=ISTEP, &
      IPRINT=IPRINT, HIST=HIST)

!                                 Compute intercept

      B0 = XMEAN(NVAR)

      DO 10  I=1, NVAR - 1

         IF (HIST(I) .GT. 0.0) THEN

            B0       = B0 - XMEAN(I)*COEF(I,1)

            J        = J + 1

            XMEAN(J) = XMEAN(I)

         END IF

   10 CONTINUE

!                                 Compute standard error of intercept

      CALL RSUBM (COVS, HIST, NIND, COVS)

      SEB0 = 1.0/NOBS + BLINF(COVS, XMEAN, XMEAN, NRA=NIND, NCA=NIND)

      SEB0 = SQRT(AOV(8)*SEB0)

!                                 Print intercept and standard error

      CALL UMACH (2, NOUT)

      WRITE (NOUT,99999) ' '

      WRITE (NOUT,99999) 'Intercept ', B0

      WRITE (NOUT,99999) 'Std. Error', SEB0

99999 FORMAT (1X, A, F10.3)

!

      END

Output

 

FORWARD SELECTION

Dependent  R-squared   Adjusted  Est. Std. Dev.
Variable   (percent)  R-squared  of Model Error
       5      97.247     96.697           2.734

                 * * * Analysis of Variance * * *
                         Sum of        Mean             Prob. of
Source           DF     Squares      Square  Overall F  Larger F
Regression        2      2641.0      1320.5    176.636    0.0000
Error            10        74.8         7.5
Total            12      2715.8

                * * * Inference on Coefficients * * *
                (Conditional on the Selected Model)
               Coef.    Standard                Prob. of    Variance
Variable    Estimate       Error  t-statistic   Larger t   Inflation
       1       1.440      0.1384       10.403     0.0000        1.06
       4      -0.614      0.0486      -12.622     0.0000        1.06

        * * * Statistics for Variables Not in the Model * * *
               Coef.    Standard  t-statistic   Prob. of    Variance
Variable    Estimate       Error     to enter   Larger t   Inflation
       2       0.416      0.1856        2.242     0.0517        18.7
       3      -0.410      0.1992       -2.058     0.0697        3.46

* * * Forward Selection Summary * * *
        Variable    Step Entered
               1             2
               4             1

Intercept    103.097
Std. Error     2.124



http://www.vni.com/
PHONE: 713.784.3131
FAX:713.781.9260