[R-meta] Non-linear model with metafor
Mustafa Kilickap
mkilickap at yahoo.com
Sun Sep 24 19:28:23 CEST 2017
Dear Wolfgang,
Thank you very much for answering to my e-mail.
I am going to explain the steps for prevalence and blood pressure data separately:
1) PREVALENCE DATA
study n nT year
a 1804 4910 2003
b 3411 10444 2007
c 8321 26499 2010
d 1647 5437 2012
n: number of hypertensive patients
nT: sample size for each study
STEPS:
library(Matrix)
library(metafor)
a<-as.data.frame(read.csv2(file.choose()))
a
min(a$year)
##min year is 2003; time is centered to 2002
a$t<-a$year-2002
a$ln.t<-log(a$t)
######EFFECT SIZES
## for raw proportion
datpr<-escalc(xi=a$n, ni=a$nT, measure="PR")
datpr
## for log transformed proportion
datpln<-escalc(xi=a$n, ni=a$nT, measure="PLN")
datpln
## for Freeman-Tukey (FT) transformed proportion
datpft<-escalc(xi=a$n, ni=a$nT, measure="PFT")
datpft
###### META-ANALYSIS
meta.pft<-rma(yi, vi, data=datpft, measure="PFT",
method="REML")
summary(meta.pft)
## I-squared 96.97%, residual heterogeneity is significant
## In order to explain the heterogeneity I planned to do a meta-regression
###### META-REGRESSION
reg<-rma(yi, vi, data=datpr, measure="PR",
method="REML", mods=~a$t)
summary(reg)
# R-squared 97.63%, I-squared 45.38%, p value for residual heterogeneity 0.166 .
# I thought that non-linear model might be more appropriate as prevalence
# of hypertension may decrease asymptotically.
# so, I tried a meta-regression for log(prevalence) on log(time) using the
#measure="PLN".
regpln_lnt<-rma(yi, vi, data=datpln, measure="PLN",
method="REML", mods=~a$ln.t)
summary(regpln_lnt)
#R-squared 100%, I-squared=0, p value for residual heterogeneity 0.66, tau=0
# beta for intercept -0.9979# beta for ln.t -0.0781
# so non-linear relationship seems to be appropriate.
# As logarithmic relationship might be difficult to grasp for readers, I
# wanted to show this relationship in another way: y= beta0* t^(beta1).
# Beta0 was obtained taking the exp. of beta for intercept:
exp(- 0.9979)
#exp(-0.9979)=0.369
#the function: 0.369*t^( - 0.0781)
# in order to show the non-linear relationship in another way, I did a
# meta-regession using the function that I obtained above.
fnx<-rma(yi, vi, data=datpr, measure="PR",
method="REML",
mods=~I(0.369*a$t^(-0.0781)))summary(fnx)
#R-squared, I-squared, and pvalue are the same above as expected.
new.t<-1:12
wi<-1/sqrt(datpr$vi)
size<- 1.5 + 3.0 * (wi - min(wi))/(max(wi) - min(wi))
plot(a$t, datpr$yi, ylim=c(0.1,0.4), xlim=c(0,12),
ylab="Prevalence",
xlab="Years (2002+)", cex=size)
pred<-predict(fnx, newmods=0.369*new.t^(-0.0781))
lines(new.t, pred$pred, col="red")
lines(new.t, pred$ci.lb, lty="dashed")
lines(new.t, pred$ci.ub, lty="dashed")
mtext(line=-3, cex=0.9, expression(italic(y==0.369(Year - 2002)
^italic(-0.0781))))
text(a$t, datpr$yi, paste(a$study), cex=0.8, pos=1)
###
MY FIRST QUESTION: Is this approach acceptable?
###
2) SYSTOLIC BLOOD PRESSURE (sbp) DATA
study n mean sd year
w 4910 130.13 22.37 2003
x 4264 127.9 19.9 2004
v 10638 126.09 16.46 2007
y 26499 120.25 24.92 2010
z 4309 133.5 27.2 2003
za 5437 127.24 18.49 2012
Assuming what I did above is true, I tried to do same thing for blood
pressure data (I mean meta-regression for log (sbp) on log (t), in order to
get the beta values). The problem is that; for proportion data “measure =PLN”
gives log-transformed prevalence and its variance. However, there is no such function
for quantitative data (only “measure = MN” is available). Even if
log-transformed blood pressure can easily be obtained, I cannot obtain the
variance of log-transformed blood pressure, so I could not do the
meta-regression.
My questions:
a) How can I get the beta values for non-linear relationship?
b) Is it acceptable to do simply lm(log(sbp)~log(t)) to get the beta
values , then test the non-linear function using these betas in metaphor as below:
a$t<-a$year-2002
a$ln.t<-log(a$t)
a$lnsbp<-log(a$mean)
kk<-lm(a$lnsbp~a$ln.t)
summary(kk)
#beta for intercept= 4.8774 ; exp(4.8793)=131.29
#beta1: -0.0266
sb<-escalc(mi=a$mean, sdi=a$sd, ni=a$n, measure="MN")
sb.reg<-rma(yi, vi, data=sb, measure="MN", method="REML", mods=~I(131.29*a$t^(-0.0266)))
c) The third question is for both prevalence and sbp data: First I got linear meta-regression (sbp on ~t; didn't show above), then non-linear meta-regression (sbp on
131.29*a$t^(-0.0266)). How can I compare these in terms of statistical significance? Is comparison based on AIC values appropriate?
Thank you very much again for your help
Kind regards,
Mustafa
On Sunday, September 24, 2017 5:33 PM, Viechtbauer Wolfgang (SP) <wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
Dear Mustafa,
I am having a hard time understanding what exactly it is that you did. Would you be able to provide a fully reproducible example (you can also use some mock data) to show your steps?
Best,
Wolfgang
-----Original Message-----
From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-project.org] On Behalf Of Mustafa Kilickap
Sent: Saturday, 23 September, 2017 11:46
To: r-sig-meta-analysis at r-project.org
Subject: [R-meta] Non-linear model with metafor
Dear metafor users,
I am doing a meta-analysis of epidemiological studies assessing prevalence of hypertension and mean systolic blood pressure (sbp). The data showed that both the prevalence and mean blood pressure decreased over 10 years. The linear models of meta-regression were significant, but residual heterogeneities were also significant. I wanted to assess a non-linear model (y=beta0*t^ (- beta1)) with metafor.
For the prevalence data, I did a meta-regression using the “measure =”PLN” on “log(year)”. The model showed I-squared =0, R-squared =100%. Then, in order to show the non-linear relationship in an understandable way, first I defined the non-linear formula by using the beta parameters estimated in this model, then I did a meta-regression on “year” (measure=PR), using the non-linear formula (mods).
For sbp data, I tried to do the same thing; in order to find the beta0 and beta1 parameters, I intended to model log(sbp) on log(t) with metafor. However, I couldn’t do the analysis as escalc function requires standard deviation (SD) of log transformed sbp, which is difficult to estimate.
I would like to know whether 1) what I did for prevalence data is correct; 2) for sbp data, whether it might be acceptable to do lm(log(sbp)~log(t)) to obtain beta parameters, then test the non-linear function in metaphor; or 3) which method is recommended.
Any help is appreciated
Kind regards
Mustafa
More information about the R-sig-meta-analysis
mailing list