tbl_summary report the z statistic for the Wilcoxon two-sample test

Following my previous question, consider:

df <- Rdatasets::rddata("military_tbl_df", package="usdatasets")
df$wage <- rnorm(nrow(df), mean = 50000, sd = 10000)

df |>
  tbl_summary(by = gender, include = c(wage, grade)) |>
  add_p() |> 
  separate_p_footnotes() |> 
  modify_column_unhide("statistic") |> 
  remove_footnote_header(columns = "statistic")

which generates

In my field is common to report the test statistic together with the p-value.

This is also suggested by Wikipedia, which provides the following example to report the result:

"Median latencies in groups E and C were 153 and 247 ms; the distributions in the two groups differed significantly (Mann–Whitney U = 10.5, n1 = n2 = 8, P < 0.05 two-tailed)."

The problem in this case is that the test statistics is very large, 122,844,026,068.

How can I have gtsummary print other type of the test statistic, like the common language effect size, the ρ statistic, or the z statistic?

According to ?tests (see: Custom Functions) you can write your own function for a test not available in gtsummary.

library(gtsummary)
packageVersion("gtsummary")
[1] ‘2.3.0’

Wilcox_z <- function(data, variable, by, ...) {
  data <- data[c(variable, by)] %>% dplyr::filter(complete.cases(.))
  
  by_x <- unique(data[[by]])[1]
  by_y <- unique(data[[by]])[2]

  x <- data[[variable]][data[[by]]==by_x]
  y <- data[[variable]][data[[by]]==by_y]
  
  res <- wilcox.test(x=x, y=y) %>%
    broom::tidy()

  library(rcompanion)
  Z <- wilcoxonZ(x=x, y=y) 
  
  res$statistic <- Z
  res
}

df |>
  tbl_summary(by = gender, include = c(wage, grade)) |>
  add_p(test = wage ~ "Wilcox_z")  |> 
  separate_p_footnotes() |> 
  modify_column_unhide("statistic") |> 
  remove_footnote_header(columns = "statistic")