[R-meta] How to interpret sigma(model)^2 in metafor

Yuhang Hu yh342 @end|ng |rom n@u@edu
Wed Jan 25 17:34:45 CET 2023


Thank you, Wolfgang. I guess my expectation was that `rstandard(res)$se^2`
should be equal to `dat$vi` since the model takes `dat$vi` as given (known
and fixed) for the sampling distribution of each residual in each row (e_ij
~ Normal[0, dat$vi])?

Thank you very much,
Yuhang

On Wed, Jan 25, 2023 at 1:40 AM Viechtbauer, Wolfgang (NP) <
wolfgang.viechtbauer using maastrichtuniversity.nl> wrote:

> Depends on what you mean by 'variance of the residuals'. With:
>
> rstandard(res)$se^2
>
> you can obtain the *sampling variance* of the residuals. Note that each
> residual has its own sampling variance.
>
> If you just want the *sample variance* of the residuals, then
>
> var(resid(res))
>
> would give you that.
>
> Best,
> Wolfgang
>
> >Thank you so much Wolfgang, for your response.
> >
> >On the same note, is there a way to extract the variance of the
> >residuals for the model, say from a fitted model like below?
> >
> >dat <- dat.konstantopoulos2011
> >res <- rma.mv(yi ~ I(year-mean(year)), vi, random = ~ 1 |district/study,
> >data= dat)
> >
> >Thank you.
> >Yuhang
> >
> >On Tue, Jan 24, 2023 at 1:21 AM Viechtbauer, Wolfgang (NP) <
> >wolfgang.viechtbauer using maastrichtuniversity.nl> wrote:
> >
> >> Dear Yuhang,
> >>
> >> sigma() is a generic function:
> >>
> >> > sigma
> >>
> >> function (object, ...)
> >>
> >> UseMethod("sigma")
> >>
> >> <bytecode: 0x2e1a3c8>
> >>
> >> <environment: namespace:stats>
> >>
> >> So, when calling sigma() on an object from the metafor package, the 'S3
> >> dispatch mechanism' will first check if there is a method for the type
> >> (i.e., class) of object that you are passing to the sigma() function.
> >>
> >> > library(metafor)
> >> > methods(sigma)
> >> [1] sigma.default* sigma.gls*     sigma.lmList*  sigma.lme*
> >> sigma.mlm*
> >> see '?methods' for accessing help and source code
> >>
> >> Since there is none (there is no sigma.rma() function or anything like
> >> it), it will call sigma.default(). So let's look at what that does:
> >>
> >> > sigma.default
> >> Error: object 'sigma.default' not found
> >>
> >> Hmmm, why can't we look at the code for this function? Note the * after
> >> sigma.default -- this indicates that the method definition is not
> >> exported.
> >> But we can still look at this with:
> >>
> >> > getAnywhere(sigma.default)
> >> A single object matching 'sigma.default' was found
> >> It was found in the following places
> >>   registered S3 method for sigma from namespace stats
> >>   namespace:stats
> >> with value
> >>
> >> function (object, use.fallback = TRUE, ...)
> >> sqrt(deviance(object, ...)/(nobs(object, use.fallback = use.fallback) -
> >>     sum(!is.na(coef(object)))))
> >> <bytecode: 0x9806a08>
> >> <environment: namespace:stats>
> >>
> >> (or, if you would happen to know that this function comes from the
> >> stats
> >> package, you could use stats:::sigma.default).
> >>
> >> So, we can see what is happening. In essence:
> >>
> >> dat <- escalc(measure="RR", ai=tpos, bi=tneg, ci=cpos, di=cneg,
> >> data=dat.bcg)
> >> res <- rma(yi, vi, data=dat)
> >> sigma(res)
> >> sqrt(deviance(res)/(nobs(res) - sum(!is.na(coef(res)))))
> >>
> >> This has, as far as I am concerned, no logical meaning for rma objects.
> >>
> >> Note that I try to be very explicit in the documentation what kind of
> >> methods are available (and meaningful) for a given object in the
> >> metafor
> >> package:
> >>
> >> https://wviechtb.github.io/metafor/reference/rma.uni.html#methods
> >>
> >> sigma() is not listed there. I cannot prevent default methods from
> >> being
> >> called unless I would actually put a sigma.rma() method into the
> >> metafor
> >> package. I have actually considered this, but I don't have a good idea
> >> what
> >> meaningful result this should return.
> >>
> >> In any case, I hope this provides you with some idea how you can dig
> >> into
> >> the code (and the mechanisms of how it is being called) in general.
> >>
> >> Best,
> >> Wolfgang
> >>
> >> >Hello Colleagues,
> >> >
> >> >By habit, I always check the variance of the residuals of my ordinary
> >> >regression models using: sigma(model)^2, which is also printed in the
> >> >output.
> >> >
> >> >I know that the variance of the residuals in meta-regression is not
> >> >estimated but rather taken as being known and fixed by virtue of user's
> >> >supplying the 'vi' or 'V' to functions such as rma.uni() and
> >> >rma.mv().
> >> >
> >> >So, I was wondering what is the interpretation of
> sigma(rma.uni_model)^2
> >> >and sigma(rma.mv_model)^2 and how they connect to the user-supplied
> >> >'vi'.
> >> >
> >> >Thank you very much for your time.
> >> >
> >> >Best,
> >> >Yuhang
>


-- 
Yuhang Hu (She/Her/Hers)
Ph.D. Student in Applied Linguistics
Department of English
Northern Arizona University

	[[alternative HTML version deleted]]



More information about the R-sig-meta-analysis mailing list