Home > Back-end >  Convert part of a statistical function's output into a data.frame
Convert part of a statistical function's output into a data.frame

Time:10-06

I was wondering if there might be a way to turn the following part of the OUTPUT of the res and res2 objects into a data.frame?

Note: answer below works with res but not res2.

A functional answer is appreciated as the data below is just toy.

library(metafor)

dat <- dat.konstantopoulos2011
res <- rma.mv(yi, vi, random = ~ 1 | district/school, data=dat)

#== OUTPUT (CAN WE TURN ONLY BELOW PART INTO A data.frame?):

#Variance Components:

#            estim    sqrt  nlvls  fixed           factor 
#sigma^2.1  0.0651  0.2551     11     no         district 
#sigma^2.2  0.0327  0.1809     56     no  district/school 

#Test for Heterogeneity:
#Q(df = 55) = 578.8640, p-val < .0001

# AND

res2 <- rma.mv(yi, vi, random = ~ factor(school) | district, data=dat)

#== OUTPUT (CAN WE TURN ONLY BELOW PART INTO A data.frame?):
#Variance Components:

#outer factor: district       (nlvls = 11)
#inner factor: factor(school) (nlvls = 11)

#            estim    sqrt  fixed 
#tau^2      0.0978  0.3127     no 
#rho        0.6653             no 

#Test for Heterogeneity:
#Q(df = 55) = 578.8640, p-val < .0001

CodePudding user response:

If there is no default/standard way to extract the data then you can manipulate the output using capture.output.

return_data <- function(res) {
  tmp <- capture.output(res)
  #data start from second line after "Variance Components:"
  start <- which(tmp == "Variance Components:")   2
  index <- which(tmp == "") 
  #Data ends before the empty line after "Variance Components:"
  end <- index[which.max(index > start)] - 1
  data <- read.table(text = paste0(tmp[start:end], collapse = '\n'), header = T)
  heterogeneity_index <- which(tmp == "Test for Heterogeneity:")   1
  list(data = data, heterogeneity = tmp[heterogeneity_index])
}

res <- rma.mv(yi, vi, random = ~ 1 | district/school, data=dat)
return_data(res)

#$data
#           estim   sqrt nlvls fixed          factor
#sigma^2.1 0.0651 0.2551    11    no        district
#sigma^2.2 0.0327 0.1809    56    no district/school

#$heterogeneity
#[1] "Q(df = 55) = 578.8640, p-val < .0001"

CodePudding user response:

Would this suit your purposes? The 'Test for Heterogeneity' doesn't really fit in the dataframe, so I added it as a seperate column and it gets duplicated as a result. I'm not sure how else you could do it.

library(tidyverse)
#install.packages("metafor")
library(metafor)
#> Loading required package: Matrix
#> 
#> Attaching package: 'Matrix'
#> The following objects are masked from 'package:tidyr':
#> 
#>     expand, pack, unpack
#> 
#> Loading the 'metafor' package (version 3.0-2). For an
#> introduction to the package please type: help(metafor)

dat <- dat.konstantopoulos2011
res <- rma.mv(yi, vi, random = ~ 1 | district/school, data=dat)
res
#> 
#> Multivariate Meta-Analysis Model (k = 56; method: REML)
#> 
#> Variance Components:
#> 
#>             estim    sqrt  nlvls  fixed           factor 
#> sigma^2.1  0.0651  0.2551     11     no         district 
#> sigma^2.2  0.0327  0.1809     56     no  district/school 
#> 
#> Test for Heterogeneity:
#> Q(df = 55) = 578.8640, p-val < .0001
#> 
#> Model Results:
#> 
#> estimate      se    zval    pval   ci.lb   ci.ub 
#>   0.1847  0.0846  2.1845  0.0289  0.0190  0.3504  * 
#> 
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

vc <- cbind(estim = res$sigma2,
            sqrt = res$sigma,
            nlvls = res$s.nlevels, 
            fixed = ifelse(res$vc.fix$sigma2, "yes", "no"), 
            factor = res$s.names,
            R = ifelse(res$Rfix, "yes", "no"),
            Test_for_heterogeneity = paste0("Q(df = ", res$k - res$p, ") = ", metafor:::.fcf(res$QE, res$digits[["test"]]), ", p-val ", metafor:::.pval(res$QEp, 
                                                                                                                               res$digits[["pval"]], showeq = TRUE, sep = " "))
            
)
rownames(vc) <- c("sigma^2.1", "sigma^2.2")
result <- as.data.frame(vc)
result
#>      estim                nlvls fixed factor            R    Test_for_heterogeneity                
#> sigma^2.1 "0.0650619442753117" "11"  "no"  "district"        "no" "Q(df = 55) = 578.8640, p-val < .0001"
#> sigma^2.2 "0.0327365170279351" "56"  "no"  "district/school" "no" "Q(df = 55) = 578.8640, p-val < .0001"

Created on 2021-10-06 by the reprex package (v2.0.1)

  • Related