Neste tutorial, vamos aprender como aplicar o método dos mínimos quadrados para fazer uma regressão linear em seus dados experimentais.
O que é Regressão Linear
Como o nome pode sugerir, é uma técnica onde tentaremos transformar uma série de dados em uma linha.
Ela é muito útil para tentarmos estimarmos ou prevermos novos valores, sem precisarmos realizar novamente o experimento.
Por exemplo, vamos supor que os pontos azuis abaixo sejam referente aos experimentos que você realizou em laboratório:
No exemplo que fizemos, usamos uma mola e várias massas. Para cada massa m, temos um tamanho de mola l. Colocamos a massa no eixo X e o comprimento no eixo Y.
Note na imagem acima, que a maioria dos pontos ficarem mais ou menos alinhados. O que a regressão linear faz é tentar achar uma reta 'ótima', mais próxima desses pontos.
Assim, você pode estimar qual seria o valor de uma mola para uma determinada massa, sem precisar fazer o experimento. Obviamente, é uma aproximação.
Método dos mínimos quadrados
Existem vários métodos e maneiras de se fazer regressão, inclusive de qualquer tipo de curva que você imaginar, não somente de retas.
O método mais conhecido, é o do mínimos quadrados. Veja que existe uma distância de cada bolinha (experimento que você fez) até nossa reta otimizada por regressão linear:
O nosso objetivo é pegar essa distância, que é o erro e minimizar o máximo possível eles, e gerar a reta de equação:
- y = A + B.x
Você pode achar a demonstração desse método em livros universitários que falam sobre Teoria dos Erros. As fórmulas para A e B são:
Dados coletados e fórmula
Para nosso tutorial, vamos gerar o gráfico para a mola menor. Queremos provar que a lei de Hooke para o sistema massa-mola é obedecido.
Equilíbrio do sistema massa-mola:
- −k∆x + mg = 0
Nos experimentos, medimos o tamanho l final da mola, que é obtido pela expressão acima:
- L = Lo + (k/g)*m
Vamos fazer a regressão linear dessa fórmula, através do método dos mínimos quadrados, com a reta:
- y = A + Bx
y ⇔ l
x ⇔ m
Regressão Linear: Fazendo os Cálculos e gráficos
Obviamente, existe tudo isso (e muito mais) prontinho em Python, para ser usado. Bastando só inserir seus dados e ele vai fazer tudo. Mas, provavelmente seu professor vai pedir que você faça isso 'na mão', pelo menos a primeira vez.
E isso é o correto. Implementar, na unha, as coisas vai te fazer um cientista ou engenheiro bem mais capacitado, vai te fazer entender melhor como as coisas funcionam e caso precise mudar algo, vai poder mexer no código e na matemática.
É muito importante, principalmente para TCC, mestrado, doutorado etc, entender bem o que está ocorrendo, pois você pode ter que fazer algo novo, para seu trabalho.
Então vamos, lá.
Primeiro, importamos a biblioteca numpy (para os cálculos) e a matplotlib.pyplot (para gerar o gráfico).
Depois, criamos dois arrays, um com as massas (massa,em gramas) e outro com o valor das distensões da mola (delta_x, em centímetros), medidos no laboratório.
Seja L0 o comprimento inicial da mola, então os comprimentos finais serão dados por L0 + delta_x
Como fizemos 7 medições, N=7. E agora é só fazer os cálculos das fórmulas de A e B.
Fizemos cada uma separadamente, para melhor entender, e depois jogamos tudo para ter as variáveis A e B. Por fim, geramos o array y com os novos dos comprimentos estimados, parar gerar a reta.
Colocamos, para aqueles que desejarem, o cálculo das incertezas de A, B, y e da constante da mola. Lembrando que o parâmetro A fica em centímetros e B fica com cm/gramas.
Como a gravidade usada é 9.8 m/s², não se esqueça de passar as outras unidades para o SI.
Nossos dados experimentais, os pontos, foram plotados com a scatter (dispersão), depois usamos uma reta pontilhada para plotar a linha da regressão, e por fim um gráfico de erros.
Note que, ao final do script, coloquei vários prints, para ir exibindo os valores das variáveis e cálculos realizados, uma etapa importante para você ir acompanhando e conferindo se está tudo ok. É muito fácil e comum cometermos erros ao longo do desenvolvimento de scripts para assuntos científicos.
Nosso script fica:
import numpy as np import matplotlib.pyplot as plt # Massas(em gramas) e respectivas Distensões da mola deltax (em centímetros) massa = np.array([10.7, 20.6, 30.8, 41.1, 50.9, 60.8, 71.2]) delta_x = np.array([6.4, 8.0, 8.8, 10.0, 11.0, 12.8, 13.5]) #Comprimentos totais da mola L0 = 5.0 L = L0 + delta_x # Seja: l = l0 + Bm a equação da reta que relaciona o tamanho # final da mola (l), com seu tamanho inicial (l0) e o valor da massa (m) # B representa a constante g/k (onde g é a acelaração da gravidade) e k # a constante da mola ## Regressão linear: y = A + Bx # x é a massa e y é o l # Número de experimentos N = 7 # Somatório de x[i] (massas) S_m = massa.sum() # Somatório de l[i] (comprimentos) S_l = L.sum() # Somatório dos quadrados x²[i] das massas S_m2 = sum(i*i for i in massa) # Somatório dos quadrados dos comprimentos l²[i] S_l2 = sum(i*i for i in L) # Somatório dos produtos de cada massa pelo seu respectivo comprimento final m[i]*l[i] S_ml = 0.0 for i in range(N): S_ml = S_ml + massa[i]*L[i] # Colocando tudo na fórmula delta = N*S_m2 - (S_m*S_m) A = ((S_m2 * S_l) - S_m * S_ml)/(delta) B = ((N*S_ml) - S_m * S_l)/(delta) # Reta com os novos valores, Yi y = np.array(A + B*massa) # Incerteza de y, A e B, com uma unidade de incerteza residuo=0.0 y=y.round(1) A=round(A,1) B=round(B,1) for i in range(N): residuo = residuo + np.power(y[i] - A - B*massa[i], 2) dy = round(np.sqrt(residuo/(N-2)),1) dA = round(dy * np.sqrt(S_m2/delta),1) dB = round(dy * np.sqrt(N/delta),1) # B e seu erro dB estão em centimetros/gramas = 10 metros/kilogramas B_kgm = B*10.0 dB_kgm = dB*10.0 # Constante da mola, e seu erro propagado de B # pois B = g/k K = 9.8/B_kgm dK = dB_kgm * ( 9.8/(B_kgm*B_kgm) ) # Gerando o gráfico fig, ax = plt.subplots() fig.suptitle('Massa(g) x Tamanho da mola (cm)') plt.rcParams.update({'font.size': 12}) ax.title.set_text('Curso Python Progressivo') ax.scatter(massa, L, marker='.', s=1) ax.plot(massa, y, color='red', linestyle='--', label='reta mmq') ax.errorbar(massa, L, yerr=dy, fmt='.k', label='dados experimentais') ax.legend(loc="lower right") plt.show() print("Somatorio das massas: ", S_m) print("Somatorio quadrado das massas: ", S_m2) print("Somatorio comprimento: ", S_l) print("Somatorio quadrado do scomprimentos: ", S_l2) print("Somatorio massa pelo comprimento menores: ", S_ml) print("Delta: ", delta) print("A: %.2f cm"% A) print("B: %.2f cm/g"% B) print("Reta da mola: y = %.1f + %.1f*m" % (A, B)) print("Erro y: ", dy) print("Erro A: ", dA) print("Erro B: ", dB) print("Constante: ", K) print("Erro de k: ", dK) print("K = %.1f +- %.1f" %(K, dK))
Gráfico da regressão linear:
Código no github:
Regressão Linear usando a biblioteca SciPy
Agora vamos deixar o Python fazer (quase) tudo por nós. Para isso, importamos a biblioteca SciPy, e desde o conjunto de funções stats, onde iremos uar a função linregress para fazer a regressão linear de nossos dados.
Lembre-se: você vai precisar de dois conjuntos de valores de mesmo tamanho, para passar para a linregress. Essa função retorna um conjunto de parâmetros, como o A e B (intercept e slope), previamente calculados, além de outros, como erros padrões.
Nosso código fica:
import numpy as np import matplotlib.pyplot as plt from scipy import stats def func(x): return result.slope * x + result.intercept # Massas(g) e respectivas Distensões da mola (cm) massa = np.array([10.7, 20.6, 30.8, 41.1, 50.9, 60.8, 71.2]) delta_x = np.array([6.4, 8.0, 8.8, 10.0, 11.0, 12.8, 13.5]) # Comprimento total da mola L0=5.0 L = np.array(L0 + delta_x) # Armazenando os parametros # Retorna a tupla com: B(coefieciente angular da reta), intercept (A, onde toca o eixo x), # rvalue (coeficiente Pearson de correlação), pvalue (valor-p da probabilidade), # stderr (erro padrão do parâmetro slope), intercept_stderr (erro padrão do parâmetro intercept) result = stats.linregress(massa, L) modelo = list(map(func, massa)) # Gerando o gráfico fig, ax = plt.subplots() fig.suptitle('Massa(g) x Tamanho da mola (cm), com SciPy') plt.rcParams.update({'font.size': 12}) ax.title.set_text('Curso Python Progressivo') ax.scatter(massa, L, marker='.') ax.plot(massa, modelo, color='red', linestyle='--', label='reta mmq') plt.show() ## Imprimindo os dados print("B (coeficiente angular): ",result.slope) print("A (onde incercepta o eixo x): ", result.intercept) print("rvalue : ", result.rvalue) print("pvalue : ", result.pvalue) print("stderr (erro padrão do slope, B): ", result.stderr) print("intercept_stderr (erro padrão do incercept, A): ",result.intercept_stderr)
E o gráfico:
Ao final, imprimimos os valores calculados pela linregress(). Veja que os valores de A e B realmente bateram com aqueles que calculamos antes. Porém, não deu exato, pois em nossos cálculos tivemos o cuidado de arredondar para uma casa decimal, pois nossas medições tinham apenas uma casa decimal. Sempre tomem cuidado com os algoritmos significativos.
