R Script for Meta-analysis of TMH
Clístenes Crístian de Carvalho
4/13/2021
R Markdown
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
SR.table<-data.frame(TP=c(11,43,36,22,41,24,34,19), FN=c(4,12,37,4,2,8,3,4), FP=c(19,30,158,3,7,9,198,2), TN=c(86,465,378,287,89,304,216,289))
SR.table$names<-c("Mostafa 2020","Panjiar 2019", "Yabuki 2019", "Rao 2018", "Nurullah 2018", "Jain 2017", "Selvi 2016", "Etezadi 2013")
SR.table$design<-c("Cohort", "Cohort", "Cohort", "Cohort", "Cohort", "Cohort", "Cohort", "Cohort")
library(mada)## Loading required package: mvtnorm## Loading required package: ellipse##
## Attaching package: 'ellipse'## The following object is masked from 'package:graphics':
##
## pairs## Loading required package: mvmeta## This is mvmeta 1.0.3. For an overview type: help('mvmeta-package').madad(SR.table)## Descriptive summary of SR.table with 8 primary studies.
## Confidence level for all calculations set to 95 %
## Using a continuity correction of 0.5 if applicable
##
## Diagnostic accuracies
## sens 2.5% 97.5% spec 2.5% 97.5%
## Mostafa 2020 0.733 0.480 0.891 0.819 0.735 0.881
## Panjiar 2019 0.782 0.656 0.871 0.939 0.915 0.957
## Yabuki 2019 0.493 0.382 0.605 0.705 0.665 0.742
## Rao 2018 0.846 0.665 0.938 0.990 0.970 0.996
## Nurullah 2018 0.953 0.845 0.987 0.927 0.857 0.964
## Jain 2017 0.750 0.579 0.867 0.971 0.946 0.985
## Selvi 2016 0.919 0.787 0.972 0.522 0.474 0.569
## Etezadi 2013 0.826 0.629 0.930 0.993 0.975 0.998
##
## Test for equality of sensitivities:
## X-squared = 43.832, df = 7, p-value = 2.3e-07
## Test for equality of specificities:
## X-squared = 543.6962, df = 7, p-value = <2e-16
##
##
## Diagnostic OR and likelihood ratios
## DOR 2.5% 97.5% posLR 2.5% 97.5% negLR 2.5% 97.5%
## Mostafa 2020 12.447 3.575 43.340 4.053 2.437 6.740 0.326 0.140 0.757
## Panjiar 2019 55.542 26.533 116.268 12.900 8.876 18.748 0.232 0.141 0.383
## Yabuki 2019 2.328 1.419 3.819 1.673 1.281 2.185 0.719 0.569 0.907
## Rao 2018 526.167 110.735 2500.120 81.795 26.223 255.137 0.155 0.063 0.383
## Nurullah 2018 260.643 51.869 1309.729 13.076 6.388 26.766 0.050 0.013 0.194
## Jain 2017 101.333 35.849 286.433 26.083 13.291 51.189 0.257 0.141 0.469
## Selvi 2016 12.364 3.738 40.889 1.921 1.672 2.208 0.155 0.052 0.462
## Etezadi 2013 686.375 118.122 3988.343 120.196 29.824 484.407 0.175 0.072 0.427
##
## Correlation of sensitivities and false positive rates:
## rho 2.5 % 97.5 %
## -0.135 -0.767 0.629par(mfrow=c(1,3))
tabletext<-cbind(SR.table$names,SR.table$design)
mada::forest(madad(SR.table),type ="sens",snames=SR.table$names)
mada::forest(madad(SR.table),type = "spec",snames=SR.table$names)
mada::forest(madad(SR.table),type = "DOR",snames=SR.table$names)## Warning in arrows(lb[i], (N:1)[i], ub[i], angle = 90, code = 3, length = 0.05, :
## zero-length arrow is of indeterminate angle and so skipped
par(mfrow=c(1,1))
mada::ROCellipse(SR.table)
#Forest Plots by univariate analysis with meta
library(meta)## Loading 'meta' package (version 4.12-0).
## Type 'help(meta)' for a brief overview.##
## Attaching package: 'meta'## The following object is masked from 'package:mada':
##
## forest
sensitivity_logit <- meta::metaprop(SR.table$TP, SR.table$TP + SR.table$FN, comb.fixed=FALSE, comb.random=TRUE, sm= "PLOGIT", method.ci= "CP", studlab=SR.table$names)
print(sensitivity_logit, digits=3)## proportion 95%-CI
## Mostafa 2020 0.733 [0.449; 0.922]
## Panjiar 2019 0.782 [0.650; 0.882]
## Yabuki 2019 0.493 [0.374; 0.613]
## Rao 2018 0.846 [0.651; 0.956]
## Nurullah 2018 0.953 [0.842; 0.994]
## Jain 2017 0.750 [0.566; 0.885]
## Selvi 2016 0.919 [0.781; 0.983]
## Etezadi 2013 0.826 [0.612; 0.950]
##
## Number of studies combined: k = 8
##
## proportion 95%-CI
## Random effects model 0.812 [0.694; 0.891]
##
## Quantifying heterogeneity:
## tau^2 = 0.6075; tau = 0.7794; I^2 = 76.5%; H = 2.06
##
## Test of heterogeneity:
## Q d.f. p-value Test
## 35.89 7 < 0.0001 Wald-type
## 44.58 7 < 0.0001 Likelihood-Ratio
##
## Details on meta-analytical method:
## - Random intercept logistic regression model
## - Maximum-likelihood estimator for tau^2
## - Logit transformation
## - Clopper-Pearson confidence interval for individual studiesforest(sensitivity_logit, digits=3, rightcols=c("effect", "ci"), xlab="Sensitivity")
specificity_logit <- meta::metaprop(SR.table$TN, SR.table$TN+ SR.table$FP, comb.fixed=FALSE, comb.random=TRUE, sm="PLOGIT", method.ci="CP", studlab=SR.table$names)
print(specificity_logit, digits=3)## proportion 95%-CI
## Mostafa 2020 0.819 [0.732; 0.887]
## Panjiar 2019 0.939 [0.915; 0.959]
## Yabuki 2019 0.705 [0.665; 0.744]
## Rao 2018 0.990 [0.970; 0.998]
## Nurullah 2018 0.927 [0.856; 0.970]
## Jain 2017 0.971 [0.946; 0.987]
## Selvi 2016 0.522 [0.472; 0.571]
## Etezadi 2013 0.993 [0.975; 0.999]
##
## Number of studies combined: k = 8
##
## proportion 95%-CI
## Random effects model 0.930 [0.809; 0.976]
##
## Quantifying heterogeneity:
## tau^2 = 2.5447; tau = 1.5952; I^2 = 98.6%; H = 8.41
##
## Test of heterogeneity:
## Q d.f. p-value Test
## 317.32 7 < 0.0001 Wald-type
## 559.54 7 < 0.0001 Likelihood-Ratio
##
## Details on meta-analytical method:
## - Random intercept logistic regression model
## - Maximum-likelihood estimator for tau^2
## - Logit transformation
## - Clopper-Pearson confidence interval for individual studiesforest(specificity_logit, digits=3, rightcols=c("effect", "ci"), xlab="Specificity")
DOR_model <- meta::metabin(TP,TP+FP,FN,FN+TN, sm="OR", comb.fixed=FALSE,comb.random=TRUE, method= "Inverse",studlab = SR.table$names, data=SR.table)
print(DOR_model)## OR 95%-CI %W(random)
## Mostafa 2020 12.4474 [ 3.5749; 43.3403] 12.5
## Panjiar 2019 55.5417 [ 26.5326; 116.2676] 13.2
## Yabuki 2019 2.3277 [ 1.4190; 3.8186] 13.4
## Rao 2018 526.1667 [110.7352; 2500.1202] 12.0
## Nurullah 2018 260.6429 [ 51.8693; 1309.7288] 11.9
## Jain 2017 101.3333 [ 35.8493; 286.4333] 12.8
## Selvi 2016 12.3636 [ 3.7384; 40.8890] 12.6
## Etezadi 2013 686.3750 [118.1219; 3988.3426] 11.6
##
## Number of studies combined: k = 8
##
## OR 95%-CI z p-value
## Random effects model 56.6070 [11.9548; 268.0386] 5.09 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 4.6250 [1.5914; 19.6191]; tau = 2.1506 [1.2615; 4.4293];
## I^2 = 94.6% [91.5%; 96.6%]; H = 4.31 [3.43; 5.41]
##
## Test of heterogeneity:
## Q d.f. p-value
## 129.82 7 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
## - Jackson method for confidence interval of tau^2 and taumeta::forest(DOR_model, digits=3, rightcols=c("effect", "ci"), xlab = "Diagnostic Odds Ratio",layout = "JAMA",leftcols = NULL)
#DSL for random effects (hetoregeneity) and MH for fixed effects (homogeneity)
(fit.DOR.DSL<-mada::madauni(SR.table,type="DOR",method = "DSL"))## Call:
## mada::madauni(x = SR.table, type = "DOR", method = "DSL")
##
## DOR tau^2
## 56.607 4.625(fit.PLR.DSL<-mada::madauni(SR.table,type="posLR"))## Call:
## mada::madauni(x = SR.table, type = "posLR")
##
## posLR tau^2
## 10.764 1.252(fit.NLR.DSL<-mada::madauni(SR.table,type="negLR"))## Call:
## mada::madauni(x = SR.table, type = "negLR")
##
## negLR tau^2
## 0.226 0.602summary(fit.DOR.DSL)## Call:
## mada::madauni(x = SR.table, type = "DOR", method = "DSL")
##
## Estimates:
## DSL estimate 2.5 % 97.5 %
## DOR 56.607 11.955 268.039
## lnDOR 4.036 2.481 5.591
## tau^2 4.625 0.000 11.912
## tau 2.151 0.000 3.451
##
## Cochran's Q: 5.692 (7 df, p = 0.576)
## Higgins' I^2: 0%summary(fit.PLR.DSL)## Call:
## mada::madauni(x = SR.table, type = "posLR")
##
## Estimates:
## DSL estimate 2.5 % 97.5 %
## posLR 10.764 4.758 24.349
## lnposLR 2.376 1.560 3.192
## tau^2 1.252 0.000 9.257
## tau 1.119 0.000 3.043
##
## Cochran's Q: 12.335 (7 df, p = 0.09)
## Higgins' I^2: 43.25%summary(fit.NLR.DSL)## Call:
## mada::madauni(x = SR.table, type = "negLR")
##
## Estimates:
## DSL estimate 2.5 % 97.5 %
## negLR 0.226 0.123 0.416
## lnnegLR -1.488 -2.099 -0.878
## tau^2 0.602 0.000 1.546
## tau 0.776 0.000 1.244
##
## Cochran's Q: 4.87 (7 df, p = 0.676)
## Higgins' I^2: 0%par(mfrow=c(1,3))
mada::forest(fit.DOR.DSL, log=FALSE)## Warning in arrows(lb[i], (N:1)[i], ub[i], angle = 90, code = 3, length = 0.05, :
## zero-length arrow is of indeterminate angle and so skippedmada::forest(fit.PLR.DSL, log=FALSE)## Warning in arrows(lb[i], (N:1)[i], ub[i], angle = 90, code = 3, length = 0.05, :
## zero-length arrow is of indeterminate angle and so skipped
## Warning in arrows(lb[i], (N:1)[i], ub[i], angle = 90, code = 3, length = 0.05, :
## zero-length arrow is of indeterminate angle and so skippedmada::forest(fit.NLR.DSL, log=FALSE)
par(mfrow=c(1,1))
#diagnostic accuracy taking into account homogeneity and heterogeneity of studies
(fit.phm.hetero<-mada::phm(SR.table,hetero=TRUE))## Call:
## phm.default(data = SR.table, hetero = TRUE)
##
## Coefficients:
## theta taus_sq
## 0.11891318 0.01580926summary(fit.phm.hetero)## Call:
## phm.default(data = SR.table, hetero = TRUE)
##
## Estimate 2.5 % 97.5 %
## theta 0.11891318 -0.11179175 0.34961811
## taus_sq 0.01580926 -0.01513495 0.04675346
##
## Log-likelihood: 10.862 on 2 degrees of freedom
## AIC: -17.7
## BIC: -17.6
##
## Chi-square goodness of fit test (Adjusted Profile Maximum Likelihood
## under heterogeneity)
##
## data: x
## Chi-square = 10.1, df = 2, p-value = 0.1205
##
##
## AUC 2.5 % 97.5 % pAUC 2.5 % 97.5 %
## 0.894 1.126 0.741 0.823 1.212 0.579(fit.phm.homo<-mada::phm(SR.table,hetero=FALSE))## Call:
## phm.default(data = SR.table, hetero = FALSE)
##
## Coefficients:
## theta
## 0.0483328summary(fit.phm.homo)## Call:
## phm.default(data = SR.table, hetero = FALSE)
##
## Estimate 2.5 % 97.5 %
## theta 0.0483328 0.03643233 0.06023326
##
## Log-likelihood: -0.847 on 1 degrees of freedom
## AIC: 3.7
## BIC: 3.8
##
## Chi-square goodness of fit test (Adjusted Profile Maximum Likelihood
## under homogeneity)
##
## data: x
## Chi-square = 41.491, df = 1, p-value = 6.516e-07
##
##
## AUC 2.5 % 97.5 % pAUC 2.5 % 97.5 %
## 0.954 0.965 0.943 0.923 0.941 0.905par(mfrow=c(1,1))
plot(fit.phm.hetero)
ROCellipse(SR.table,add = TRUE)
(fit.reitsma<-mada::reitsma(SR.table))## Call: reitsma.default(data = SR.table)
##
## Fixed-effects coefficients:
## tsens tfpr
## (Intercept) 1.3833 -2.5265
##
## 8 studies, 2 fixed and 3 random-effects parameters
## logLik AIC BIC
## 14.6673 -19.3346 -15.4717summary(fit.reitsma) ## Call: reitsma.default(data = SR.table)
##
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: REML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) 1.383 0.324 4.272 0.000 0.749 2.018 ***
## tfpr.(Intercept) -2.526 0.599 -4.218 0.000 -3.701 -1.352 ***
## sensitivity 0.800 - - - 0.679 0.883
## false pos. rate 0.074 - - - 0.024 0.205
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.771 1.000 .
## tfpr 1.651 -0.202 1.000
##
## logLik AIC BIC
## 14.667 -19.335 -15.472
##
## AUC: 0.902
## Partial AUC (restricted to observed FPRs and normalized): 0.855
##
## HSROC parameters
## Theta Lambda beta sigma2theta sigma2alpha
## 0.149 3.751 0.762 0.508 3.059plot(fit.reitsma,sroclwd = 2)
points(fpr(SR.table),sens(SR.table), pch=2)
legend("bottomright", c("data","summary estimate"),pch=c(2,1))
library(metafor)## Loading required package: Matrix## Loading 'metafor' package (version 2.4-0). For an overview
## and introduction to the package please type: help(metafor).##
## Attaching package: 'metafor'## The following objects are masked from 'package:meta':
##
## baujat, forest, funnel, funnel.default, labbe, radial, trimfill## The following object is masked from 'package:mada':
##
## forest## The following object is masked from 'package:mvmeta':
##
## bluplibrary(metaviz)
library(dmetar)## Extensive documentation for the dmetar package can be found at:
## www.bookdown.org/MathiasHarrer/Doing_Meta_Analysis_in_R/legend("bottomleft",c("SROC", "conf. region"), lwd=c(2,1))
dat<-metafor::escalc(measure = "OR", ai=TP,bi= FP,ci=FN,di=TN,data = SR.table)
funnel_TMH<-metafor::rma.uni(yi,vi,data=dat,method="REML")
metaviz::funnelinf(funnel_TMH,n=1,xlab = "DOR",ylab = "Standard error of DOR", egger=TRUE, trim_and_fill = TRUE) 
dmetar::eggers.test(x=DOR_model)## Warning in dmetar::eggers.test(x = DOR_model): Your meta-analysis contains k =
## 8 studies. Egger's test may lack the statistical power to detect bias when the
## number of studies is small (i.e., k<10).## Intercept ConfidenceInterval t p
## Egger's test 7.654 3.146-12.162 3.379 0.01488#summary sensitivity and specificity for studies with common threshold
SR.table.2<-data.frame(TP=c(43,3,22,41,24,34,19), FN=c(12,3,4,2,8,3,4), FP=c(30,191,3,7,9,198,2), TN=c(465,412,287,89,304,216,289))
(fit.reitsma.2<-mada::reitsma(SR.table.2))## Call: reitsma.default(data = SR.table.2)
##
## Fixed-effects coefficients:
## tsens tfpr
## (Intercept) 1.5575 -2.6653
##
## 7 studies, 2 fixed and 3 random-effects parameters
## logLik AIC BIC
## 14.5936 -19.1872 -15.9919summary(fit.reitsma.2) ## Call: reitsma.default(data = SR.table.2)
##
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: REML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) 1.557 0.261 5.958 0.000 1.045 2.070 ***
## tfpr.(Intercept) -2.665 0.683 -3.904 0.000 -4.003 -1.327 ***
## sensitivity 0.826 - - - 0.740 0.888
## false pos. rate 0.065 - - - 0.018 0.210
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.437 1.000 .
## tfpr 1.762 0.180 1.000
##
## logLik AIC BIC
## 14.594 -19.187 -15.992
##
## AUC: 0.887
## Partial AUC (restricted to observed FPRs and normalized): 0.856
##
## HSROC parameters
## Theta Lambda beta sigma2theta sigma2alpha
## 0.900 4.455 1.395 0.454 1.263#summary sensitivity and specificity for studies with prespecified threshold
SR.table.3<-data.frame(TP=c(3,22,41,24,34), FN=c(3,4,2,8,3), FP=c(191,3,7,9,198), TN=c(412,287,89,304,216))
(fit.reitsma.3<-mada::reitsma(SR.table.3))## Call: reitsma.default(data = SR.table.3)
##
## Fixed-effects coefficients:
## tsens tfpr
## (Intercept) 1.6654 -2.2436
##
## 5 studies, 2 fixed and 3 random-effects parameters
## logLik AIC BIC
## 8.7411 -7.4822 -5.9693summary(fit.reitsma.3) ## Call: reitsma.default(data = SR.table.3)
##
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: REML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) 1.665 0.459 3.627 0.000 0.766 2.565 ***
## tfpr.(Intercept) -2.244 0.830 -2.704 0.007 -3.870 -0.617 **
## sensitivity 0.841 - - - 0.683 0.929
## false pos. rate 0.096 - - - 0.020 0.350
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.824 1.000 .
## tfpr 1.823 0.030 1.000
##
## logLik AIC BIC
## 8.741 -7.482 -5.969
##
## AUC: 0.912
## Partial AUC (restricted to observed FPRs and normalized): 0.871
##
## HSROC parameters
## Theta Lambda beta sigma2theta sigma2alpha
## 0.484 3.986 0.794 0.773 2.912#sensitivity analyses
#external laryngeal manipulation
SR.table$burp<-c("yes","yes", "yes", "no", "yes", "yes", "no", "yes")
SR.table$patient<-c("no","yes", "yes", "yes", "yes", "yes", "yes", "no")
fit.burp <- mada::reitsma(SR.table,formula = cbind(tsens, tfpr) ~ burp,method = "ml")
summary(fit.burp)## Call: reitsma.default(data = SR.table, formula = cbind(tsens, tfpr) ~
## burp, method = "ml")
##
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: ML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) 2.029 0.571 3.556 0.000 0.911 3.147 ***
## tsens.burpyes -0.903 0.643 -1.404 0.160 -2.163 0.357
## tfpr.(Intercept) -2.149 1.104 -1.947 0.052 -4.313 0.015 .
## tfpr.burpyes -0.483 1.273 -0.379 0.704 -2.978 2.012
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.569 1.000 .
## tfpr 1.509 -0.397 1.000
##
## logLik AIC BIC
## 15.613 -17.227 -11.819#risk of bias for patient selection
fit.patient <- mada::reitsma(SR.table,formula = cbind(tsens, tfpr) ~ patient,method = "ml")
summary(fit.patient)## Call: reitsma.default(data = SR.table, formula = cbind(tsens, tfpr) ~
## patient, method = "ml")
##
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: ML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) 1.287 0.641 2.009 0.045 0.031 2.543 *
## tsens.patientyes 0.099 0.728 0.136 0.892 -1.328 1.526
## tfpr.(Intercept) -3.096 1.127 -2.747 0.006 -5.304 -0.887 **
## tfpr.patientyes 0.763 1.291 0.591 0.554 -1.767 3.293
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.707 1.000 .
## tfpr 1.508 -0.241 1.000
##
## logLik AIC BIC
## 14.757 -15.514 -10.106#index test
SR.table$index<-c("yes","yes", "no", "no", "yes", "no", "no", "yes")
fit.index <- mada::reitsma(SR.table,formula = cbind(tsens, tfpr) ~ index,method = "ml")
summary(fit.index)## Call: reitsma.default(data = SR.table, formula = cbind(tsens, tfpr) ~
## index, method = "ml")
##
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: ML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) 1.126 0.394 2.858 0.004 0.354 1.899 **
## tsens.indexyes 0.481 0.577 0.834 0.404 -0.650 1.613
## tfpr.(Intercept) -2.171 0.764 -2.840 0.005 -3.670 -0.673 **
## tfpr.indexyes -0.689 1.089 -0.633 0.527 -2.823 1.445
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.652 1.000 .
## tfpr 1.493 -0.169 1.000
##
## logLik AIC BIC
## 15.006 -16.012 -10.603#reference standardbaujat_TMH<-metafor::baujat(res,xlim=c(0,58),ylim = c(0,43),symbol = "slab",ylab="Influence on Overall Result")
SR.table$standard<-c("unclear","yes", "no", "yes", "unclear", "unclear", "unclear", "yes")
fit.standard <- mada::reitsma(SR.table,formula = cbind(tsens, tfpr) ~ standard,method = "ml")
summary(fit.standard)## Call: reitsma.default(data = SR.table, formula = cbind(tsens, tfpr) ~
## standard, method = "ml")
##
## Bivariate diagnostic random-effects meta-regression
## Estimation method: ML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept) -0.027 0.264 -0.104 0.917 -0.546 0.491
## tsens.standardunclear 1.685 0.384 4.388 0.000 0.932 2.437 ***
## tsens.standardyes 1.402 0.371 3.781 0.000 0.675 2.128 ***
## tfpr.(Intercept) -0.872 1.050 -0.831 0.406 -2.930 1.186
## tfpr.standardunclear -0.990 1.182 -0.838 0.402 -3.307 1.326
## tfpr.standardyes -3.049 1.247 -2.446 0.014 -5.492 -0.606 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Variance components: between-studies Std. Dev and correlation matrix
## Std. Dev tsens tfpr
## tsens 0.123 1.000 .
## tfpr 1.046 1.000 1.000
##
## logLik AIC BIC
## 22.150 -26.299 -19.346#Baujat Plot
par(mar=c(5,4,2,2))
res<-metafor::rma(yi,vi,data=dat,method = "FE",slab = names)
baujat_TMH<-metafor::baujat(res,xlim=c(0,58),ylim = c(0,43),symbol = "slab",ylab="Influence on Overall Result")
baujat_TMH## x y ids slab
## Mostafa 2020 0.09642097 0.007155780 1 Mostafa 2020
## Panjiar 2019 11.85819206 2.909609521 2 Panjiar 2019
## Yabuki 2019 55.07824454 43.080173368 3 Yabuki 2019
## Rao 2018 19.89229953 0.921460945 4 Rao 2018
## Nurullah 2018 11.92110650 0.512985826 5 Nurullah 2018
## Jain 2017 12.83380494 1.419501454 6 Jain 2017
## Selvi 2016 0.11218485 0.009117155 7 Selvi 2016
## Etezadi 2013 18.02938605 0.648598567 8 Etezadi 2013Including Plots
You can also embed plots, for example:
Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.