Gotmadstacks - 10 months ago 72

R 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)`

`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

`b:c`

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`

Answer Source

`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
```