I compiled 11 correlations from k=11 studies. Let's just call the two variables $A$ and $B$. The correlations were very high and ranged in between $r=.54-.84$, but there is one study that reported a perfect correlation of $r=1.00$. My problem is when trying to aggregate effect sizes in a random effects model (weighting by inverse variance) the correlation of $r=1.00$ can not be used since any correlation of anything greater than or equal to 1 will have 0 variance. I have also tried the following:

  • I z-transformed but the estimate is now infinity.
  • I set the correlation to be $r=.9999$ but then the aggregated correlation sky rockets from $\bar{r}=.73$ to $\bar{r}=.98$. Also the amount of decimal places I add greatly impacts the aggregate correlation.

I do not want to remove the correlation from my meta-analysis as it fits the inclusion criteria. But I am not sure how to add it in to the meta-analysis in a way that makes sense


1 Answer 1


Leaving aside whether a correlation of 1 is really possible or not (since I know nothing about the context here), let me suggest some possible approaches. I'll simulate some data (in R) for illustration purposes:

ri <- c(runif(8, .54, .84), .54, .84, 1)
ni <- rep(20, 11)

And I will use the metafor package for the following computations:


First, we can just treat the correlations as they are and compute the corresponding sampling variances in the usual manner:

dat <- escalc(measure="COR", ri=ri, ni=ni)

This yields:

       yi     vi
1  0.5741 0.0237
2  0.7267 0.0117
3  0.7228 0.0120
4  0.7270 0.0117
5  0.7983 0.0069
6  0.7321 0.0113
7  0.5428 0.0262
8  0.6098 0.0208
9  0.5400 0.0264
10 0.8400 0.0046
11 1.0000 0.0000

Since the sampling variance of a correlation coefficient is usually computed based on the equation $$Var[r_i] = v_i = \frac{(1-r_i^2)^2}{n_i-1}$$ you will get a variance of 0 for the study where $r = 1$. But we can still go ahead and meta-analyze these studies with a standard random-effects model (using REML estimation by default):

rma(yi, vi, data=dat)

This yields:

Random-Effects Model (k = 11; tau^2 estimator: REML)

tau^2 (estimated amount of total heterogeneity): 0.0156 (SE = 0.0118)
tau (square root of estimated tau^2 value):      0.1250

Model Results:

estimate      se     zval    pval   ci.lb   ci.ub
  0.7524  0.0497  15.1449  <.0001  0.6550  0.8498  ***

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

A few notes about this:

  1. You will not get the results from the $Q$-test for heterogeneity since it cannot be computed when at least one sampling variance is equal to 0.
  2. Depending on the actual data, it may happen that $\tau^2$ (the amount of heterogeneity in the true correlation coefficients) is estimated to be 0. In this case, model fitting will fail, since the weights in a random-effects model are $w_i = 1 / (v_i + \hat{\tau}^2)$ and if $\hat{\tau}^2 = 0$ then you get division by zero. One could then switch to another estimator for $\tau^2$ that avoids this problem. For example, the 'Sidik-Jonkman' estimator is strictly positive and you can use this one with rma(yi, vi, data=dat, method="SJ").
  3. If $\hat{\tau}^2$ is close to zero, then a lot of weight will be placed on the study with $r = 1$. This can drag up the pooled estimate quite a bit, which may be undesirable. A possible solution is to choose custom weights that will avoid this problem. For example, with rma(yi, vi, data=dat, method="SJ", weights=ni) we use the strictly positive SJ estimator and use custom weights $w_i = n_i$, so studies are weighted by their sample sizes only, which may counteract this issue (assuming that the study with $r = 1$ doesn't also have a relative large sample size compared to the rest of the studies).

Usually, we would want to apply the r-to-z transformation when meta-analyzing correlation coefficients. However, if you try:

dat <- escalc(measure="ZCOR", ri=ri, ni=ni)

you will find (as noted in the question) that this leads to a $z_r$ value of infinity for the study with $r = 1$. One could indeed use some kind of ad-hoc adjustment:

ri[ri == 1] <- .9999
dat <- escalc(measure="ZCOR", ri=ri, ni=ni)

This yields:

       yi     vi
1  0.6536 0.0588
2  0.9217 0.0588
3  0.9134 0.0588
4  0.9224 0.0588
5  1.0938 0.0588
6  0.9332 0.0588
7  0.6082 0.0588
8  0.7085 0.0588
9  0.6042 0.0588
10 1.2212 0.0588
11 4.9517 0.0588

Since the sampling variance of an r-to-z transformed correlation coefficient is estimated with $$Var[z_{r_i}] = v_i = \frac{1}{n_i-3}$$ all variances all positive (and they all happen to be the same here, since I simulated the data this way, but this is not relevant). Now we can proceed with the analysis:

res <- rma(yi, vi, data=dat)

This yields:

Random-Effects Model (k = 11; tau^2 estimator: REML)

tau^2 (estimated amount of total heterogeneity): 1.5042 (SE = 0.6990)
tau (square root of estimated tau^2 value):      1.2265
I^2 (total heterogeneity / total variability):   96.24%
H^2 (total variability / sampling variability):  26.57

Test for Heterogeneity:
Q(df = 10) = 265.7141, p-val < .0001

Model Results:

estimate      se    zval    pval   ci.lb   ci.ub
  1.2302  0.3770  3.2635  0.0011  0.4914  1.9690  **

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The pooled estimate still needs to be back-transformed, which we can do with:

predict(res, transf=transf.ztor)

This yields:

  pred  ci.lb  ci.ub   pi.lb  pi.ub
0.8426 0.4553 0.9618 -0.8577 0.9989

With this approach, the choice of the adjustment (do you change $r = 1$ to .9999 or .99 or something else?) can have quite a large influence on the results.

I cannot tell you which of these approaches is "correct". You could also use and report multiple approaches as sensitivity analyses. Most importantly, you should be transparent about what you did.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.