Coefficient table does not have NA rows in rank-deficient fit; how to insert them?

Question
library(lmPerm)
x <- lmp(formula = a ~ b * c + d + e, data = df, perm = "Prob")

summary(x)  # truncated output, I can see `NA` rows here!

#Coefficients: (1 not defined because of singularities)
#                 Estimate Iter Pr(Prob)
#b                   5.874   51    1.000
#c                -30.060   281    0.263
#b:c                   NA    NA       NA
#d1               -31.333    60    0.633
#d2                33.297   165    0.382
#d3               -19.096    51    1.000
#e                  1.976    NA       NA

I want to pull out the Pr(Prob) results for everything, but

y <- summary(x)$coef[, "Pr(Prob)"]

#(Intercept)           b            c           d1           d2 
# 0.09459459  1.00000000   0.26334520   0.63333333   0.38181818 
#         d3           e 
# 1.00000000          NA 

This is not what I want. I need b:c row, too, in the right position.

An example of the output I would like from the above would be:

# (Intercept)           b            c    b:c           d1           d2 
#  0.09459459  1.00000000   0.26334520     NA   0.63333333   0.38181818 
#         d3            e 
# 1.00000000           NA 

I also would like to pull out the Iter column that corresponds to each variable. Thanks.


Show source
| R   | permutation   | regression   | linear-regression   | lm   2016-11-28 13:11 1 Answers

Answers ( 1 )

  1. 2016-11-28 14:11

    lmp is based on lm and summary.lmp also behaves like summary.lm, so I will first use lm for illustration, then show that we can do the same for lmp.


    lm and summary.lm

    Have a read on ?summary.lm and watch out for the following returned values:

    coefficients: a p x 4 matrix with columns for the estimated
                  coefficient, its standard error, t-statistic and
                  corresponding (two-sided) p-value.  Aliased coefficients are
                  omitted.
    
         aliased: named logical vector showing if the original coefficients are
                  aliased.
    

    When you have rank-deficient models, NA coefficients are omitted in the coefficient table, and they are called aliased variables. Consider the following small, reproducible example:

    set.seed(0)
    zz <- xx <- rnorm(10)
    yy <- rnorm(10)
    fit <- lm(yy ~ xx + zz)
    
    coef(fit)  ## we can see `NA` here
    #(Intercept)          xx          zz 
    #  0.1295147   0.2706560          NA 
    
    a <- summary(fit)  ## it is also printed to screen
    #Coefficients: (1 not defined because of singularities)
    #            Estimate Std. Error t value Pr(>|t|)
    #(Intercept)   0.1295     0.3143   0.412    0.691
    #xx            0.2707     0.2669   1.014    0.340
    #zz                NA         NA      NA       NA
    
    b <- coef(a)  ## but no `NA` returned in the matrix / table
    #             Estimate Std. Error   t value  Pr(>|t|)
    #(Intercept) 0.1295147  0.3142758 0.4121051 0.6910837
    #xx          0.2706560  0.2669118 1.0140279 0.3402525
    
    d <- a$aliased
    #(Intercept)          xx          zz 
    #      FALSE       FALSE        TRUE 
    

    If you want to pad NA rows to coefficient table / matrix, we can do

    ## an augmented matrix of `NA`
    e <- matrix(nrow = length(d), ncol = ncol(b),
                dimnames = list(names(d), dimnames(b)[[2]]))
    ## fill rows for non-aliased variables
    e[!d] <- b
    
    #             Estimate Std. Error   t value  Pr(>|t|)
    #(Intercept) 0.1295147  0.3142758 0.4121051 0.6910837
    #xx          0.2706560  0.2669118 1.0140279 0.3402525
    #zz                 NA         NA        NA        NA
    

    lmp and summary.lmp

    Nothing needs be changed.

    library(lmPerm)
    fit <- lmp(yy ~ xx + zz, perm = "Prob")
    a <- summary(fit)  ## `summary.lmp`
    b <- coef(a)
    
    #              Estimate Iter  Pr(Prob)
    #(Intercept) -0.0264354  241 0.2946058
    #xx           0.2706560  241 0.2946058
    
    d <- a$aliased
    #(Intercept)          xx          zz 
    #      FALSE       FALSE        TRUE 
    
    e <- matrix(nrow = length(d), ncol = ncol(b),
                dimnames = list(names(d), dimnames(b)[[2]]))
    e[!d] <- b
    
    #              Estimate Iter  Pr(Prob)
    #(Intercept) -0.0264354  241 0.2946058
    #xx           0.2706560  241 0.2946058
    #zz                  NA   NA        NA
    

    If you, want to extract Iter and Pr(Prob), just do

    e[, 2]  ## e[, "Iter"]
    #(Intercept)          xx          zz 
    #        241         241          NA 
    
    e[, 3]  ## e[, "Pr(Prob)"]
    #(Intercept)          xx          zz 
    #  0.2946058   0.2946058          NA 
    
◀ Go back