My question is a bit more general, and was already asked here enter link description here. However, there seems to be no general solution. Therefore, I try to use an example. What is an efficient way to proceed to get the summary function applied to a glht object print the degrees of freedom used by the glht to compute the p-values? I think I should add an argument to the summary function, but which one?
The most reliable way to find the possible arguments is to look into the function code. The summary function returns a list object, what governs what's printed is the print()
function, so that's what we want to analize.
But print()
is a generic function, this is an Object Oriented Programming that basically means "a function that has different methods for each type of object it can receive". Running the code print
gets us:
function (x, ...)
UseMethod("print")
So we can't access the specific method of print()
that it's used with a glht summary. To do that, we can either use debugging, or use the sloop
package:
sloop::s3_methods_generic('print')
# A tibble: 403 × 4
generic class visible source
<chr> <chr> <lgl> <chr>
1 print aareg FALSE registered S3method
2 print abbrev FALSE registered S3method
3 print acf FALSE registered S3method
4 print all_vars FALSE registered S3method
5 print anova FALSE registered S3method
6 print Anova FALSE registered S3method
7 print anova.loglm FALSE registered S3method
8 print any_vars FALSE registered S3method
9 print aov FALSE registered S3method
10 print aovlist FALSE registered S3method
# ℹ 393 more rows
These are all the methods that you have installed for print()
. Now select the ones that have something to do with glht:
sloop::s3_methods_generic('print') |> filter(grepl('glht', class))
# A tibble: 3 × 4
generic class visible source
<chr> <chr> <lgl> <chr>
1 print confint.glht FALSE registered S3method
2 print glht FALSE registered S3method
3 print summary.glht FALSE registered S3method
Thus, we want the code for the method print.summary.glht()
(this is the default syntax for methods, function.class):
sloop::s3_get_method('print.summary.glht')
function (x, digits = max(3, getOption("digits") - 3), ...)
{
cat("\n\t", "Simultaneous Tests for General Linear Hypotheses\n\n")
if (!is.null(x$type))
... #I omitted the rest of the function
Now we can see that the arguments are x
, digits
(that are globally defined), and ...
(which aren't used anywhere in the method). So, in conclusion, there isn't an argument to display the degrees of freedom.
But, now that you have access to the function's code, you can alter it yourself, adding cat("Degrees of Freedom: ", x$df, "\n")
near the end, and saving it with the same name:
print.summary.glht <- function (x, digits = max(3, getOption("digits") - 3), ...){
cat("\n\t", "Simultaneous Tests for General Linear Hypotheses\n\n")
if (!is.null(x$type))
cat("Multiple Comparisons of Means:", x$type, "Contrasts\n\n\n")
call <- if (isS4(x$model))
x$model@call
else x$model$call
if (!is.null(call)) {
cat("Fit: ")
print(call)
cat("\n")
}
pq <- x$test
mtests <- cbind(pq$coefficients, pq$sigma, pq$tstat, pq$pvalues)
error <- attr(pq$pvalues, "error")
pname <- switch(x$alternative, less = paste("Pr(<", ifelse(x$df ==
0, "z", "t"), ")", sep = ""), greater = paste("Pr(>",
ifelse(x$df == 0, "z", "t"), ")", sep = ""), two.sided = paste("Pr(>|",
ifelse(x$df == 0, "z", "t"), "|)", sep = ""))
colnames(mtests) <- c("Estimate", "Std. Error", ifelse(x$df ==
0, "z value", "t value"), pname)
type <- pq$type
if (!is.null(error) && error > .Machine$double.eps) {
sig <- which.min(abs(1/error - (10^(1:10))))
sig <- 1/(10^sig)
}
else {
sig <- .Machine$double.eps
}
cat("Linear Hypotheses:\n")
alt <- switch(x$alternative, two.sided = "==", less = ">=",
greater = "<=")
rownames(mtests) <- paste(rownames(mtests), alt, x$rhs)
printCoefmat(mtests, digits = digits, has.Pvalue = TRUE,
P.values = TRUE, eps.Pvalue = sig)
cat("Degrees of Freedom: ", x$df, "\n") #our new line
switch(type, univariate = cat("(Univariate p values reported)"),
`single-step` = cat("(Adjusted p values reported -- single-step method)"),
Shaffer = cat("(Adjusted p values reported -- Shaffer method)"),
Westfall = cat("(Adjusted p values reported -- Westfall method)"),
cat("(Adjusted p values reported --", type, "method)"))
cat("\n\n")
invisible(x)
}
Example of outcome:
lmod <- lm(Fertility ~ ., data = swiss)
K <- diag(length(coef(lmod)))[-1,]
rownames(K) <- names(coef(lmod))[-1]
glht(lmod, linfct = K) |> summary() |> print()
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = Fertility ~ ., data = swiss)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
Agriculture == 0 -0.17211 0.07030 -2.448 0.0793 .
Examination == 0 -0.25801 0.25388 -1.016 0.7847
Education == 0 -0.87094 0.18303 -4.758 <0.001 ***
Catholic == 0 0.10412 0.03526 2.953 0.0233 *
Infant.Mortality == 0 1.07705 0.38172 2.822 0.0325 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Degrees of Freedom: 41 #our new line
(Adjusted p values reported -- single-step method)