Correlações¶
Criamos umas sequências de dados artificiais ...
import numpy as np
# data (imaginary!)
N = 100 # sample dimension
X, Y = [], [] # observables
for n in range(N):
x = np.random.normal(0,1)# random d's between -1 and 1
y = 50 * x + 30 * np.random.normal(0,1)
X, Y = X + [x] , Y + [y]
Um scatter plot dos dados é obtido usando o método scatter da interface matplotlib.pyplot da livraria matplotlib
from matplotlib import pyplot as plt
ax = plt.scatter(X, Y, s=20, label='data')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()
O coeficiente de correlação de Pearson e o $p$-value podem ser calculados usando o método pearson do módulo scipy.stats
import scipy.stats as stat
R = stat.pearsonr(X,Y)
R
O coeficiente de correlação de Spearman e o $p$-value podem ser calculados usando o método spearmanr do módulo scipy.stats
import scipy.stats as stat
rho = stat.spearmanr(X,Y)
rho
Outra possibilidade é introduzir diretamente os dados como séries, e usar a função OLS (Ordinary Least Squares, ou seja, mínimos quadrados) do módulo statsmodels, que produz uma tabela de resultados mais completa.
Por exemplo, na segunte tabela, produzida por Galileo, foram registadas as alturas iniciais $x$ e as distâncias atingidas $y$ de um objeto lançado com a mesma força horizontal.
import statsmodels.api as sm
X = [100, 200, 300, 450, 600, 800, 1000] # valores das alturas iniciais
Y = [253, 337, 395, 451, 495, 534, 574 ] # valores das distâncias atingidas
model = sm.OLS(Y, sm.add_constant(X), missing='drop') # modelo Y=hX+b, os dados em faltas (com NANs) são ignorados
results = model.fit()
print(results.summary()) # print results inline
results.summary()
O valor de $R^2$, assim como o $p$-value, parecem sugerir uma boa correlação linear, embora um gráfico sugere claramente uma lei diferente da lei linear.
from seaborn import regplot
ax3 = regplot(x=X, y=Y)
Mais uma possibilidade é, naturalmente, importar um ficheiro de dados num DataFrame...
Por exemplo o seguinte ficheiro, de Our World in Data, relaciona a pobreza (definida como a percentagem de população que vive com menos de 3.1 dolares por dia) e a expetativa de vida.
import pandas as pd
df = pd.read_csv('/Users/sal/Lectures/mqq/Data/poverty-vs-life-expectancy.csv')
df
df2000 = df[df["Year"] > 2013 ]
df2000
import statsmodels.api as sm
L = df2000.iloc[:,5].values # life expectancy
P = df2000.iloc[:,6].values # poverty
model = sm.OLS(L, sm.add_constant(P), missing='drop') # model L=hP+b, os dados em faltas (com NANs) são ignorados
results = model.fit()
print(results.summary()) # print results inline
results.summary()
from seaborn import regplot
ax3 = regplot(x=P, y=L)
