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

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 at 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

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Wednesday, 25 January, 2023 17:35
To: Viechtbauer, Wolfgang (NP)
Cc: r-sig-meta-analysis at r-project.org
Subject: Re: Re: [R-meta] How to interpret sigma(model)^2 in metafor

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 at 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 at 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

They are not quite the same. dat$vi are not the sampling variances of the
residuals, but of the sampling errors. Those are indeed assumed to be known
and fixed. However, the sampling variances of the residuals depend on the
model you are fitting.

Best,
Wolfgang

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Wednesday, 25 January, 2023 17:35
To: Viechtbauer, Wolfgang (NP)
Cc: r-sig-meta-analysis at r-project.org
Subject: Re: Re: [R-meta] How to interpret sigma(model)^2 in metafor

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 at 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 at 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

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Wednesday, 25 January, 2023 19:13
To: Viechtbauer, Wolfgang (NP)
Cc: r-sig-meta-analysis at r-project.org
Subject: Re: Re: [R-meta] How to interpret sigma(model)^2 in metafor

Thank you,?Wolfgang. I guess there is a gap in my understanding, then, regarding
the difference?between sampling errors which I think for the current model (i.e.,
yij = B0 + B1*yearij + uj + vij + eij) are defined as:

eij = yij - theta_ij (true effect_ij), with Var(eij) = dat$vi for each eij ~
Normal[0, dat$vi].

and residuals i.e., Var(yij - fitted(yij))?.

I thought that the variances of the sampling errors (i.e., dat$vi) which indicate
that effect sizes collected from the literature are estimated with error
eventually are the same as the Var(yij - fitted(yij)) thinking that fitted(yij)
is taken as theta_ij.

Best,
Yuhang

On Wed, Jan 25, 2023 at 9:48 AM Viechtbauer, Wolfgang (NP)
<wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
They are not quite the same. dat$vi are not the sampling variances of the
residuals, but of the sampling errors. Those are indeed assumed to be known and
fixed. However, the sampling variances of the residuals depend on the model you
are fitting.

Best,
Wolfgang

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Wednesday, 25 January, 2023 17:35
To: Viechtbauer, Wolfgang (NP)
Cc: r-sig-meta-analysis at r-project.org
Subject: Re: Re: [R-meta] How to interpret sigma(model)^2 in metafor

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 at 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 at 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