# How to meta-analyze latent correlations when r=1?

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

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:

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


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

library(metafor)


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)
dat


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


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)
dat


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)
dat


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)
res


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.