Pearson Correlation in Python

Question

I have written the following code and I am new to Python. Please let me know if I can make improvements to the following program:

import math

X = [12,11,13,13, 9,10,10,13, 5,10,10,13,10,10, 5, 8, 9, 8, 8, 9, 9,10,11, 5,12]
Y = [11,10,10,10, 9,13,10,11, 6, 7,13,14,14,11,11,10,10, 7, 8,12,11,11, 8, 7,13]

Xmean = sum(X)/len(X)
Ymean = sum(Y)/len(Y)

x = [var-Xmean for var in X]
y = [var-Ymean for var in Y]

xy =[a*b for a,b in list(zip(x,y))]
sum_xy = sum(xy)

x_square = [a*a for a in x]
y_square = [b*b for b in y]

sum_x_square = sum(x_square)
sum_y_square = sum(y_square)

sum_x_square_sum_y_square = sum_x_square*sum_y_square
sqrt_sum_x_square_sum_y_square = math.sqrt(sum_x_square_sum_y_square)

r = sum_xy/sqrt_sum_x_square_sum_y_square

print(r)

The above code is written for this formula :

Answer 1

Note: I did not check whether your code matches the formula correctly and I am not an expert of the mathematical background. I only refer to your code.

Style

[12,11,13,13, 9,10,10,13, 5,10,10,13,10,10, 5, 8, 9, 8, 8, 9, 9,10,11, 5,12]

This looks very inconsistent. You should add a space after each comma.

Naming

Xmean and X are unusual names for variables. x_square is the style you should use for Python, which is lower case with underscores, also known as "snake case".

Your names however make no sense at all semantically. You can name local variables x and y if they are e. g. coordinates, or if you use them for iterating. If they represent lists of numbers, why not just call them numbers1 and numbers2? (This might still not be perfect, but it shows what the variables mean. Edit: Note that this mainly concerns the problem of having two variables that only differ in being lower/upper case, and not nowing what x is somewhere inside your code, whereas numbers tells you that it is a collection or sequence of numbers. See below.) Note that I used a lower case x to refer to your X. Having another x later is the next problem. How do you know the difference without reading your whole code again?

x = [var-Xmean for var in X]

This took me some time to understand. You can improve it by using a common iterator variable like i and adding spaces around operators. With these improvements, and also better variable names, this could be for example:

subtracted_average = [i - mean_x for i in numbers1]

sum_xy = sum(xy)

Saving intermediate results of calculations like this is good, when you need them more than once. But here you don't. So sum(xy) is more readible if it is inside of another formula than sum_xy.

sum_x_square_sum_y_square

If you read this, do you know what it means without reading the code next to it? Is it a sum of sums? A product of sums? Or a completely different thing altogether? When trying to find a good variable names, think about "what is this"? The answer might be e. g. "the product of sums of squared numbers". product_of_sums_of_squared_numbers is long, but it tells you what it is.

Implementing Formulas

Judge for yourself. What is more readible?

1.

sum_x_square = sum(x_square)
sum_y_square = sum(y_square)

sum_x_square_sum_y_square = sum_x_square*sum_y_square
sqrt_sum_x_square_sum_y_square = math.sqrt(sum_x_square_sum_y_square)

or 2.

result = math.sqrt(sum(x_square) * sum(y_square))

The second one is actually how you write formulas in maths. So when implementing a function to represent that formula, go the same way. (Hint: Maybe this code should be a function.)

Using Functions

def pearson_correlation(numbers_x, numbers_y):
    mean_x = sum(numbers_x)/len(numbers_x)
    mean_y = sum(numbers_y)/len(numbers_y)

    subtracted_mean_x = [i - mean_x for i in numbers_x]
    subtracted_mean_y = [i - mean_y for i in numbers_y]

    x_times_y = [a * b for a, b in list(zip(subtracted_mean_x, subtracted_mean_y))]

    x_squared = [i * i for i in numbers_x]
    y_squared = [i * i for i in numbers_y]

    return sum(x_times_y) / math.sqrt(sum(x_squared) * sum(y_squared))

Call it like this:

# Contrary to what I said above about variables in general,
# the names X and Y are ok here because they match the names in the formula.
X = [12, 11, 13, 13, 9, 10, 10, 13, 5, 10, 10, 13, 10, 10, 5, 8, 9, 8, 8, 9, 9, 10, 11, 5, 12]
Y = [11, 10, 10, 10, 9, 13, 10, 11, 6, 7, 13, 14, 14, 11, 11, 10, 10, 7, 8, 12, 11, 11, 8, 7, 13]

print(pearson_correlation(X, Y))

Answer 2

Since you are doing statistics, you could familiarize yourself with numpy. It excels at numerical calculations by doing them in C in the background. Also, most operations and functions are implemented in a vectorized way, which allows you to get rid of all of your explicit loops and makes it faster at the same time.

Here is the code from @D. Everhard's answer directly re-written using numpy:

import numpy as np

def pearson_correlation(numbers_x, numbers_y):
    mean_x = numbers_x.mean()
    mean_y = numbers_y.mean()
    subtracted_mean_x = numbers_x - mean_x
    subtracted_mean_y = numbers_y - mean_y
    x_times_y = subtracted_mean_x * subtracted_mean_y
    x_squared = subtracted_mean_x**2
    y_squared = subtracted_mean_y**2
    return x_times_y.sum() / np.sqrt(x_squared.sum() * y_squared.sum())

You could even make this almost a one-liner by foregoing the intermediate variables:

def pearson_correlation(numbers_x, numbers_y):
    x = numbers_x - numbers_x.mean()
    y = numbers_y - numbers_y.mean()
    return (x * y).sum() / np.sqrt((x**2).sum() * (y**2).sum())

Note that for this to work, you need to have numbers_x and numbers_y defined as numpy.array:

numbers_x = np.array([1., 2., ...])