3
$\begingroup$

I am trying to study the correlation between 2 variables non-normally distributed. These are the values (n. 134) of measured level of a substance (bisphenol A) in urine in mothers and relative newborns: bpaTm (Bisphenol A Total Mother) and bpaTn (Bisphenol A Total Newborn) (see supplementary at the bottom for the original variables). Here their distribution plot:

bpaTm distribution: enter image description here

bpaTn distribution: enter image description here

These are not normalized neither by log() nor boxcox() transformation (see supplementary data at the bottom for their distribution after transformation and boxcox() transformation meaning). But the values of Pearson correlation test changes with the transformation (I guess variables are getting more "linearly" correlated):

> cor.test( db0$bpaTm, db0$bpaTn )

    Pearson's product-moment correlation

data:  db0$bpaTm and db0$bpaTn
t = 0.077775, df = 132, p-value = 0.9381
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.1630063  0.1761554
sample estimates:
        cor 
0.006769264 

> 
> cor.test( log(db0$bpaTm), log(db0$bpaTn) )

    Pearson's product-moment correlation

data:  log(db0$bpaTm) and log(db0$bpaTn)
t = -1.5421, df = 132, p-value = 0.1254
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.29593959  0.03740469
sample estimates:
       cor 
-0.1330276 

> 
> cor.test( boxcox_transform(db0$bpaTm), boxcox_transform(db0$bpaTn) )

    Pearson's product-moment correlation

data:  boxcox_transform(db0$bpaTm) and boxcox_transform(db0$bpaTn)
t = 34.013, df = 132, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.9267169 0.9623740
sample estimates:
      cor 
0.9474099 

my question is: it is correct to apply a Pearson correlation on such transformed data? Or there are some kind of theoretical errors in this transformation?

PS: I know I can use Spearman and Kendall's tau in order to have a rank comparison not influenced by transformation, but I would like to understand the "correctness" of transformations over the Pearson correlation that is influenced by these transformations

SUPPLEMENTARY DATA:

Here the plot of the tranformed variables:

log(bpaTm) distribution: enter image description here log(bpaTn) distribution: enter image description here boxcox_transform(bpaTm) distribution: enter image description here boxcox_transform(bpaTn) distribution: enter image description here

and here the variables data:

> dput( db0$bpaTm )
c(0.0801, 0.0062, 0.0445, 0.0033, 0.1039, 0.1584, 0.0426, 0.004, 
0.1322, 0.0021, 0.0773, 1.0562, 0.0303, 0.0429, 0.6408, 0.2638, 
0.0065, 2.7172, 0.1071, 0.2451, 0.3235, 0.0547, 0.084, 0.3314, 
0.0079, 0.004, 0.0027, 0.0996, 0.0022, 0.6869, 2.3778, 0.3234, 
0.0977, 0.0798, 0.5955, 0.2373, 2.7369, 0.0179, 0.4769, 0.8089, 
0.1191, 0.4511, 0.1964, 0.0025, 0.0034, 0.002, 0.0362, 0.1084, 
0.0814, 0.0052, 0.1043, 0.0055, 8e-04, 0.0062, 0.0026, 0.0034, 
0.0044, 1.1179, 0.0029, 0.0018, 0.002, 0.7946, 0.0337, 0.0914, 
0.0793, 0.0026, 0.0037, 0.0015, 0.011, 0.0201, 0.0064, 0.077, 
0.7026, 0.708, 0.3732, 0.0092, 0.7062, 0.3412, 0.2587, 0.0049, 
0.2711, 0.4696, 0.2196, 0.5888, 0.2751, 0.0057, 0.3905, 0.7022, 
0.0049, 0.5612, 0.0014, 0.0013, 0.3752, 0.349, 0.0011, 0.0138, 
1.4973, 0.2776, 0.4582, 0.0072, 0.5075, 0.0022, 0.7452, 0.2239, 
0.8222, 0.1791, 0.4244, 0.0042, 0.007, 0.085, 0.3591, 0.9713, 
7e-04, 0.0066, 0.2405, 0.0013, 0.2325, 0.6365, 0.0027, 0.3721, 
0.0045, 0.0046, 0.1868, 0.3672, 0.0018, 0.0718, 0.0727, 0.0446, 
0.1136, 0.0283, 0.1102, 0.1274, 0.0261, 0.0844)

c(0.1273, 0.01858, 0.02705, 0.02037, 0.84932, 0.02425, 0.42594, 
0.11851, 0.23176, 0.55762, 0.022, 0.01888, 0.02659, 0.02928, 
0.05532, 0.03981, 0.27789, 2.07942, 0.01527, 0.43818, 0.08772, 
0.01087, 0.49211, 0.01918, 0.02688, 0.03996, 0.13425, 0.03981, 
0.23498, 0.14297, 1.40088, 0.01869, 0.09098, 0.66275, 0.0149, 
1.19778, 0.01913, 0.05903, 0.38875, 0.02282, 0.01389, 0.03707, 
0.03262, 0.08987, 0.01933, 0.02119, 0.14978, 0.03883, 0.30822, 
0.94456, 0.01933, 0.3146, 0.02581, 0.01678, 0.03177, 5.41625, 
0.01925, 0.20909, 0.54781, 0.26473, 0.03725, 0.05871, 0.02637, 
0.20771, 0.02161, 0.07835, 0.01947, 0.69323, 2.20418, 0.03079, 
0.83059, 0.01598, 0.02812, 0.02298, 0.02378, 0.74808, 0.05344, 
0.15396, 0.233, 0.14283, 0.0232, 0.40811, 0.16316, 0.09577, 0.52982, 
0.04868, 0.28942, 0.34124, 0.01783, 0.01883, 0.02244, 0.78543, 
0.16619, 0.01928, 0.3058, 3.39427, 0.17005, 0.02649, 0.05704, 
0.12563, 0.4699, 0.37497, 0.02692, 0.19658, 3.83661, 0.02296, 
0.58335, 0.03506, 3.00608, 0.01256, 0.03751, 0.02485, 0.73003, 
4.07519, 5.66704, 7.75299, 0.43675, 0.01065, 0.31909, 0.09174, 
0.02113, 0.0218, 0.02975, 0.0263, 1.90869, 0.01766, 0.05318, 
0.03136, 1.1696, 0.78481, 0.30664, 0.21008, 0.123, 1.34231)

boxcox() transformation is:

boxcox_transform <- function( var1 ) {
    Box = boxcox(var1 ~ 1,              # Transform as a single vector
                 lambda = seq(-6,6,0.1)      # Try values -6 to 6 by 0.1
                 )

    Cox = data.frame(Box$x, Box$y)            # Create a data frame with the results
    Cox2 = Cox[with(Cox, order(-Cox$Box.y)),] # Order the new data frame by decreasing y
#Cox2[1,]                                  # Display the lambda with the greatest log likelihood
lambda = Cox2[1, "Box.x"]                 # Extract that lambda
var1_box = (db0$bpsFn ^ lambda - 1)/lambda   # Transform the original data
    return(var1_box)
}

Here the scatter plots of the variables:

non transformed enter image description here log transformed enter image description here boxocox transformed enter image description here

$\endgroup$
3
  • $\begingroup$ nice question, can you put the scatter plots as well? $\endgroup$
    – TPArrow
    Dec 10, 2019 at 11:44
  • $\begingroup$ @TPArrow edited with scatter plots, any idea? Thank you in advance $\endgroup$
    – cccnrc
    Dec 10, 2019 at 12:31
  • $\begingroup$ I do not see any difference between BC and log transformed scatterplots. $\endgroup$
    – TPArrow
    Dec 11, 2019 at 15:26

2 Answers 2

3
$\begingroup$

To do a correct analysis on empirically transformed variables you are compelled to account for model uncertainty when assessing the final correlation coefficient. For example if you bootstrapped the whole process you'd see a lot of transformation uncertainty propagate into the confidence interval for the correlation coefficient. See for example The Cost of Data Analysis by Faraway.

Far, far easier is to use a robust correlation index that not only does not assume a distribution for either variable, but does not need transformations (if they are monotonic) nor an assumption of linearity. I'd just go with Spearman's $\rho$ here. Use your time wisely.

$\endgroup$
2
  • $\begingroup$ Thank you very much for the answer, do you have the article PDF you cited? I cannot access it and I would love to read it $\endgroup$
    – cccnrc
    Dec 10, 2019 at 14:43
  • $\begingroup$ This looks like a preprint of the final paper: people.bath.ac.uk/jjf23/papers/cda.pdf $\endgroup$ Dec 10, 2019 at 14:49
1
$\begingroup$

I try to look at the problem from a different angle.

Here the data violates two fundamental assumptions of the Pearson cor.test,

  1. Your variables should be approximately normally distributed
  2. There should be no significant outliers

See here for more details. Then cor.test result is not reliable for your data.

Coming back to the effect of transformations on the (linear) correlation. In this answer, I only explain log transformation and others could comment or improve.

Note that log is a monotone one-to-one transformation, as a result, the linear correlation should not change under this transformation. See below that we generate 10 linearly correlated observations.

set.seed(123456)
x = c(runif(10, 1, 2))
y = x * 2 + rnorm(10, sd = .1)
plot(x, y)

enter image description here

and

> cor(x,y)
[1] 0.9849312
> cor(log(x),log(y))
[1] 0.9822018

that shows no sign of significant change in the correlations. Now we add an outlier to this data,

set.seed(123456)
x = c(runif(10, 1, 2), 10)
y = c(head(x * 2 + rnorm(11, sd = .1), 10), 3)
plot(x, y)

enter image description here

The correlation coefficient of the raw and transformed data are listed below,

> cor(x, y)
[1] 0.1170776
> cor(log(x), log(y))
[1] 0.3400024

That shows some strong signals of the change in the correlation.

Conclusion: in your data, the log transformation acts toward reduces the effect of outliers (and shrinking the long tails of the distribution) and this improves the correlation coefficient.

$\endgroup$

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.