Cálculo Numérico Faculdade de Engenharia, Arquiteturas e Urbanismo – FEAU Prof. Dr. Sergio Pilling (IPD/ Física e Astronomia) V – Ajuste de curvas pelo método dos mínimos quadrados Objetivos: O objetivo desta aula é apresentar o método dos mínimos quadrados (MMQ) como outra forma de aproximação de funções. Ao contrário do polinômio interpolador visto no capitulo anterior, agora não é necessário que o ajuste passe exatamente por cima dos pontos ajustados. Em outras palavras, com esse método encontramos uma função ϕ(x) de um certo tipo pré-estabelecido (ex. reta, parábola, senoide) que melhor ajusta um conjunto de pontos ou uma função dada. 1. Introdução Como vimos na última aula, uma forma de se trabalhar com uma função definida por uma tabela de valores é a interpolação. Contudo, a interpolação pode não ser aconselhável quando: 1) É preciso obter um valor aproximado da função em algum ponto fora do intervalo de tabelamento (extrapolação). 2) Os valores tabelados são resultado de experimentos físicos, pois estes valores poderão conter erros inerentes que, em geral, não são previsíveis. Surge então a necessidade de se ajustar a estas funções tabeladas uma função que seja uma “boa aproximação” para as mesmas e que nos permita “extrapolar” com certa margem de segurança. Assim, o objetivo deste processo é aproximar uma função f(x) por outra função ϕ(x), escolhida de uma família de funções ou por uma soma de funções em duas situações distintas: Domínio discreto: quando a função f é dada por uma tabela de valores. Domínio contínuo: quando a função f é dada por sua forma analítica. Veremos nesta aula o método de ajuste de curva aos pontos experimentais (caso discreto) pelo método dos mínimos quadrados! V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 1 2 Caso Discreto O problema do ajuste de curvas no caso em que se tem uma tabela de m pontos com x1, x2, x3 , … , xm ∈[a,b], consiste em: “escolhidas” n funções contínuas g1(x), g2 (x), g3(x), … , gn(x), contínuas em [a,b], obter n constantes a1, a2, a3, …, an tais que a função ϕ(x) = a1 g1(x) + a2 g2 (x)+ a3 g3 (x)+ … + an gn (x) se aproxime ao máximo de f (x). Este modelo matemático é linear pois os coeficientes que devem ser determinados a1, a2, a3, …, an aparecem linearmente, embora as funções g1(x), g2(x), g3(x), …, gn(x) possam ser funções não lineares de x, como por exemplo, g1(x)= x2, g2(x)= ex, g3(x)= (1+x)2, etc. Surge então a primeira pergunta: Como escolher as funções contínuas g1(x), g2 (x), g3(x), … , gn (x) ? Esta escolha pode ser feita observando o gráfico dos pontos tabelados (diagrama de dispersão) ou baseando-se em fundamentos teóricos do experimento que forneceu a tabela. Portanto, dada uma tabela de pontos (x1, f(x1)), (x2,f(x2)), ...., (xn,f(xn)), deve-se, em primeiro lugar colocar estes pontes num gráfico cartesiano e a partir daí pode-se visualizar a curva que melhor se ajusta aos dados. EXEMPLO 1 Seja a tabela de pontos abaixo: O diagrama de dispersão para esses pontos é apresentado ao lado: Esse diagrama se assemelha muito a uma parábola com centro na origem, não e? Portanto, nesse caso, é natural escolhermos apenas uma função g1(x)=x2 e procurarmos então ϕ(x) = ax2 (equação geral de uma parábola passando pela origem). V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 2 EXEMPLO 2 Se considerarmos uma experiência onde foram medidos vários valores de corrente elétrica (i) que passa por uma resistência (R) submetida a várias tensões (V), colocando os valores correspondentes de corrente elétrica e tensão em um gráfico podemos ter a figura ao lado: Neste caso, existe uma fundamentação teórica relacionando a corrente com a tensão (V= R i; Lei de Ohm), isto é, V é uma função linear de i. Assim, g1(i)= i e ϕ(i)= a g1(i) = a i. Queremos ajustar nesse caso uma reta. Surge agora a segunda pergunta: Qual parábola com equação αx2 melhor se ajusta ao diagrama do exemplo 1 e qual reta, passando pela origem, melhor se ajusta ao diagrama do exemplo 2? No caso geral, escolhidas as funções g1(x), g2(x), ..., gn(x), temos de estabelecer o conceito de proximidade entre as funções ϕ(x) e f(x) para obter as constantes a1, a2, a3, …, an. Uma idéia é impor que o desvio entre f(x) e ϕ(x), ou seja, dk=(f(xk)-ϕ(xk)) seja mínimo para todos os pontos (k =1, 2, ...., m). Existem varias formas de impor que os desvios sejam mínimos. Veremos nessa aula o método dos mínimos quadrados. ϕ(x) Seja dk = f (xk) − ϕ (xk) o desvio em xk . a a a a ϕ A derivada tem que ser igual a zero! a a a a ϕ V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 3 Para isto é necessário que: aj aj a a a a a a a a Obs: A derivada tem que ser zero para acharmos o valor mínimo de F. a a a a a a a a a a a a a a a a a a a a a a a a a ^ ^ a a1 a2 a1 a2 an an aa11 a2 an ^ b V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 4 ai a1, a2, ..., an Leitura opcional: produto escalar ^ ^ ai b^i Produto escalar a1 a2 an a1 a2 a3 a1 a2 a3 an an OBS: ^ ^ ai b^i V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 5 Na prática, o funcionamento do MMQ pode ser dividido em 4 passos: PASSO 1 – Depois de escolhida a função ajuste ϕ(x) identificar nela as funções auxiliares g(x) tal que ϕ(x) seja do tipo: n φ ( x ) = ∑ ai gi ( x ) = a1 g1 ( x ) + a2 g 2 ( x ) + a3 g3 ( x )... + an g n ( x ) n∈Ι i =1 PASSO 2 – Montar o sistema de equações. O numero de equações do sistema igual ao numero de funções auxiliares gi(x) (igual ao numero de incógnitas ai) Ex 1. No caso da reta: ϕ (x) = a1 + a2 x→ g1(x) = 1 e g2(x) = x Teremos um sistema com 2 equações. ⎛ a11 ⎜⎜ ⎝ a21 a12 ⎞ ⎛ a1 ⎞ ⎛ b1 ⎞ ⎟×⎜ ⎟ = ⎜ ⎟ a22 ⎟⎠ ⎜⎝ a2 ⎟⎠ ⎜⎝ b2 ⎟⎠ incógnitas Ex 2. No caso de uma parábola: ϕ (x) = a1 + a2 x + a3x2 → g1(x) = 1 , g2(x) = x e g3(x) = x2 Teremos um sistema com 3 equações. ⎛ a11 ⎜ ⎜ a21 ⎜a ⎝ 31 a12 a22 a32 a13 ⎞ ⎛ a1 ⎞ ⎛ b1 ⎞ ⎟ ⎜ ⎟ ⎜ ⎟ a23 ⎟ × ⎜ a2 ⎟ = ⎜ b2 ⎟ a33 ⎟⎠ ⎜⎝ a3 ⎟⎠ ⎜⎝ b3 ⎟⎠ incógnitas Ex 3. No caso de uma exponencial simples: ϕ (x) = a1 ex → g1(x) = ex Teremos um sistema com 1 equação. a11 × a1 = b1 incógnita PASSO 3 – Calcular os coeficientes aij e bi do passo 2. Esses coeficientes são definidos pelos seguintes somatórios e após seu calculo obteremos números. número de pontos experimentais. m aij = ∑ g i ( xk ) g j ( xk ) =a ji k =1 m bi = ∑ f ( xk ) g i ( xk ) k =1 PASSO 4 – Reescrever o sistema de equações do passo 2 (agora os aij e bi são números) e resolvêlo, por exemplo, utilizando o método de eliminação de Gauss ou algum método iterativo (Gauss-Jacobi ou Gauss-Seidel). V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 6 Exercício 1 ϕ (x) = a1 + a2 x→ g1(x) = 1 e g2(x) = x ϕ (x) = a1 + a2 x + a3x2 → g1(x) = 1 , g2(x) = x e g3(x) = x2 Solução a) Nesse caso temos f(x) ≈ ϕ (x) = a1 + a2x o que resulta em termos g1(x) = 1 e g2(x) = x Para encontrarmos a1 e a2 resolveremos o sistema de 2 equações abaixo: ⎛ a11 ⎜⎜ ⎝ a21 a12 ⎞ ⎛ a1 ⎞ ⎛ b1 ⎞ ⎟×⎜ ⎟ = ⎜ ⎟ a22 ⎟⎠ ⎜⎝ a2 ⎟⎠ ⎜⎝ b2 ⎟⎠ Escrevendo o sistema em termos dos aij e bi , ficamos assim: Temos 8 pontos experimentais 8 ⎡8 ⎤ ⎡8 ⎤ ⎢∑ g1 ( x k ) g1 ( x k )⎥a1 + ⎢∑ g 2 ( x k ) g1 ( x k )⎥a 2 = ∑ f ( x k ) g1 ( x k ) k =1 ⎣ k =1 ⎦ ⎣ k =1 ⎦ 8 ⎡8 ⎤ ⎡8 ⎤ ⎢∑ g1 ( xk ) g 2 ( xk )⎥a1 + ⎢∑ g 2 ( xk ) g 2 ( xk )⎥a2 = ∑ f ( xk ) g 2 ( xk ) k =1 ⎣ k =1 ⎦ ⎣ k =1 ⎦ Cada somatório da parte esquerda resultará em: 8 8 ∑ g1 ( xk ) g1 ( xk ) = ∑ (g1 ( xk )) = 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 = 8 2 k =1 k =1 1 1 8 ∑ g ( x ) g ( x ) = 1 ×1 + 2 ×1 + 3 ×1 + 4 ×1 + 5 ×1 + 6 ×1 + 7 ×1 + 8 ×1 = 36 k 2 k =1 k 1 xk 1 8 ∑ g ( x ) g ( x ) = 1 ×1 + 1 × 2 + 1 × 3 + 1 × 4 + 1 × 5 + 1 × 6 + 1 × 7 + 1 × 8 = 36 k =1 k 1 k 2 xk 1 8 8 ∑ g ( x ) g ( x ) = ∑ (g ( x )) 2 k =1 k 2 xk k 2 k =1 2 k = 1 × 1 + 2 × 2 + 3 × 3 + 4 × 4 + 5 × 5 + 6 × 6 + 7 × 7 + 8 × 8 = 204 xk V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 7 Cada somatório da parte direita resultará em: 8 ∑ f (x k =1 8 ∑ f (x k =1 1 k ) g1 ( xk ) = 0.5 × 1 + 0.6 × 1 + 0.9 × 1 + 0.8 × 1 + 1.2 × 1 + 1.5 × 1 + 1.7 × 1 + 2.0 × 1 = 9.2 1 k xk ) g 2 ( xk ) = 0.5 × 1 + 0.6 × 2 + 0.9 × 3 + 0.8 × 4 + 1.2 × 5 + 1.5 × 6 + 1.7 × 7 + 2.0 × 8 = 50.5 xk Reescrevendo o sistema de equações teremos: 8a1 + 36a2 = 9.2 - 36a1 - 162a2 = -41.5 L1=-4.5L1 36a1 + 204a2 = 50.5 36a1 + 204a2 = 50.5 Subtraindo as duas equações encontramos: 50.5 − 41.5 50 − 204 × 0.214 a2 = = 0.214 e a1 = = 0.176 204 − 162 36 Podemos agora escrever a equação que ajusta os pontos experimentais f(x) ≈ϕ (x) = a1 + a2x. Resposta: ϕ(x) = 0.176 + 0.214 x Solução b) Nesse caso temos f(x) ≈ ϕ (x) = a1 + a2 x + a3 x2 o que resulta em termos g1(x) = 1, g2(x) = x e g3(x) = x2. De forma análoga ao caso anterior, para encontrarmos a1, a2 e a3 resolveremos o sistema de 3 equações abaixo: ⎛ a11 ⎜ ⎜ a21 ⎜a ⎝ 31 a12 a22 a32 a13 ⎞ ⎛ a1 ⎞ ⎛ b1 ⎞ ⎟ ⎜ ⎟ ⎜ ⎟ a23 ⎟ × ⎜ a2 ⎟ = ⎜ b2 ⎟ a33 ⎟⎠ ⎜⎝ a3 ⎟⎠ ⎜⎝ b3 ⎟⎠ Escrevendo o sistema em termos dos aij e bi , ficamos assim: 8 ⎡8 ⎤ ⎡8 ⎤ ⎡8 ⎤ ⎢∑ g1 ( xk ) g1 ( xk )⎥a1 + ⎢∑ g 2 ( xk ) g1 ( xk )⎥a2 + ⎢∑ g 3 ( xk ) g1 ( xk )⎥a3 = ∑ f ( xk ) g1 ( xk ) k =1 ⎣ k =1 ⎦ ⎣ k =1 ⎦ ⎣ k =1 ⎦ 8 8 8 8 ⎡ ⎤ ⎡ ⎤ ⎡ ⎤ ⎢∑ g1 ( xk ) g 2 ( xk )⎥a1 + ⎢∑ g 2 ( xk ) g 2 ( xk )⎥a2 + ⎢∑ g 3 ( xk ) g 2 ( xk )⎥a3 = ∑ f ( xk ) g 2 ( xk ) k =1 ⎣ k =1 ⎦ ⎣ k =1 ⎦ ⎣ k =1 ⎦ 8 ⎡8 ⎤ ⎡8 ⎤ ⎡8 ⎤ g ( x ) g ( x ) a + g ( x ) g ( x ) a + g ( x ) g ( x ) a = ⎢ ∑ 1 k 3 k ⎥ 1 ⎢ ∑ 2 k 3 k ⎥ 2 ⎢ ∑ 3 k 3 k ⎥ 3 ∑ f ( xk ) g 3 ( x k ) k =1 ⎣ k =1 ⎦ ⎣ k =1 ⎦ ⎣ k =1 ⎦ V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 8 Cada somatório da parte esquerda resultará em: 8 8 ∑ g1 ( xk ) g1 ( xk ) = ∑ (g1 ( xk )) = 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 + 1 × 1 = 8 k =1 1 8 2 k =1 1 ∑ g ( x ) g ( x ) = 1 ×1 + 2 ×1 + 3 ×1 + 4 ×1 + 5 ×1 + 6 ×1 + 7 ×1 + 8 ×1 = 36 2 k =1 k 1 xk 8 ∑ g (x 3 k =1 k 1 ) g1 ( xk ) = 12 × 1 + 2 2 × 1 + 32 × 1 + 4 2 × 1 + 52 × 1 + 6 2 × 1 + 7 2 × 1 + 82 × 1 = 204 k xk2 1 8 ∑ g ( x ) g ( x ) = 1 ×1 + 1 × 2 + 1 × 3 + 1 × 4 + 1 × 5 + 1 × 6 + 1 × 7 + 1 × 8 = 36 1 k =1 k 2 xk 1 8 k 8 ∑ g 2 ( xk ) g 2 ( xk ) = ∑ (g2 ( xk )) = 1 ×1 + 2 × 2 + 3 × 3 + 4 × 4 + 5 × 5 + 6 × 6 + 7 × 7 + 8 × 8 = 204 2 k =1 8 k =1 xk xk ∑ g ( x )g ( x ) = 1 2 3 k =1 2 k k xk2 8 × 1 + 2 2 × 2 + 32 × 3 + 4 2 × 4 + 52 × 5 + 6 2 × 6 + 7 2 × 7 + 82 × 8 = 1296 xk ∑ g ( x ) g ( x ) = 1×1 1 k =1 k 2 k xk2 1 8 3 ∑ g ( x )g ( x ) = 1×1 2 2 k =1 8 k k =1 k xk xk2 xk2 xk2 ∑ g (x 3 3 8 + 1 × 2 2 + 1 × 32 + 1 × 4 2 + 1 × 52 + 1 × 6 2 + 1 × 7 2 + 1 × 82 = 204 + 2 × 2 2 + 3 × 32 + 4 × 4 2 + 5 × 52 + 6 × 6 2 + 7 × 7 2 + 8 × 82 = 1296 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 k ) g 3 ( x k ) = ∑ (g 3 ( x k ) ) =1 × 1 + 2 × 2 + 3 × 3 + 4 × 4 + 5 × 5 + 6 × 6 + 7 × 7 + 8 × 8 = 8772 2 k =1 Cada somatório da parte direita resultara em: 8 ∑ f (x k =1 8 ∑ f (x k =1 8 ∑ k =1 k ) g1 ( xk ) = 0.5 × 1 + 0.6 × 1 + 0.9 × 1 + 0.8 × 1 + 1.2 × 1 + 1.5 × 1 + 1.7 × 1 + 2.0 × 1 = 9.2 1 k ) g 2 ( xk ) = 0.5 × 1 + 0.6 × 2 + 0.9 × 3 + 0.8 × 4 + 1.2 × 5 + 1.5 × 6 + 1.7 × 7 + 2.0 × 8 = 50.5 xk f ( x k ) g 3 ( x k ) = 0.5 × 12 + 0.6 × 2 2 + 0.9 × 32 + 0.8 × 4 2 + 1.2 × 5 2 + 1.5 × 6 2 + 1.7 × 7 2 + 2.0 × 8 2 = 319.1 xk2 Reescrevendo o sistema de equações teremos: 8a1 + 36a2 + 204a3 = 9.2 36a1 + 204a2 + 1296a3 = 50.5 204a1 +1296a2 + 8772a3 = 319.1 V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 9 Nesse caso utilizaremos o método direto de eliminação de Gauss para resolver o sistema de equações. matriz sanduíche otimizada 204 1296 8772 36 204 1296 8 36 204 1ª etapa de eliminação 319.1 50.5 9.2 2ª etapa de eliminação 204 1296 8772 319.1 0 -24.706 -252 -5.812 0 0 11.193 0.354 204 1296 8772 0 -24.706 -252 0 -14.823 -140 319.1 -5.812 -3.133 re-escrevendo o sistema de equações 204a1 + 1296 a2 + 8772a3 = 319.1 -24.706 a2 - 252 a3 = -5.812 11.193a3 = 0.354 Resolvendo o sistema de baixo para cima encontramos: a3 = 0.0316, a2 = -0.0871 e a1 = 0.7587 Podemos agora escrever a equação da parábola que melhor ajusta os pontos experimentais f(x) ≈ ϕ (x) = a1 + a2x + a3x2. Resposta: ϕ (x) = 0.7587 -0.0871 x + 0.0316 x2 Exercício 2 Resolveremos agora o exemplo 1 que vimos no inicio da aula. A partir da função tabelada abaixo, desenhamos o diagrama de dispersão e percebemos que a melhor curva que ajusta os pontos seria um a parábola passando pela origem, ou seja, f(x) ≈ ϕ(x)= a1x2 (neste caso teremos apenas 1 função g1(x) = x2). V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 10 Como só temos uma função de g(x) e f(x) ≈ ϕ(x)= a1x2 temos de resolver apenas a equação e com isso encontramos diretamente o valor de a1 a11 × a1 = b1 11 ∑ g (x 1 k =1 xk2 11 k ) g1 ( xk ) × a1 = ∑ g1 ( xk ) f ( xk ) k =1 xk2 xk2 Resolvendo os dois somatórios temos: 11 ∑ (g ( x ) ) 2 k =1 k 1 = 1 + 0.3164 + 0.1296 + 0.0625 + 0.0081 + 0 + 0.0016 + 0.0256 + 0.0625 + 0.2401 + 1 = 2.8464 (xk2)2 11 ∑ f (x k =1 k ) g1 ( xk ) = 2.05 + 0.6486 + 0.162 + 0.1 + 0.045 + 0 + 0.008 + 0.096 + 0.128 + 0.588 + 2.05 = 5.8756 xk2 Logo nossa equação é 2.8464 a = 5.8756 a = 2.0642 Então ϕ(x) = 2.0642 x é a parábola que melhor se aproxima dos pontos tabelados segundo o método dos mínimos quadrados. 2 2 ϕ(x) = 2.0642 x a>2.0642 a<2.0642 V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 11 Exercício proposto 1 ϕ (x) = a1 + a2 x→ g1(x) = 1 e g2(x) = x ϕ (x) = a1 + a2 x + a3x2 → g1(x) = 1 , g2(x) = x e g3(x) = x2 x 0 1 2 3 4 y 27 42 60 87 127 Resp: a) ϕ(x) = 19,6 + 24,5 x b) ϕ(x) = 28,02 + 7,64 x + 4,21 x2 Exercício proposto 2 DICA: Para usar o método dos mínimos quadrados é necessário termos ϕ(x) no formato abaixo: n f ( x ) ≈ φ ( x ) = ∑ ai g i ( x ) = a1 g1 ( x ) + a2 g 2 ( x ) + a3 g 3 ( x )... + an g n ( x ) n∈Ι i =1 Portanto temos que reescrever a equação proposta para o ajuste y = a ebx Aplicando ln dos dois lados temos: ln(y) = ln(a ebx) = ln(a) + ln (ebx) = ln(a) + bx Fazendo ln(y) = y* e ln(a) = a* ficamos com a equação da reta ao lado: y* = a* + bx Basta agora reescrever a tabela acima usando x e y* = ln(y) e aplicar o método MMQ. Depois de encontrarmos os valores a*= ln(a) e b escrevemos a função original y = a ebx Resp: a) y = 3,469 + 0,355 x → y= 32,104 e * 0,355 x ln ⎛⎜ y a ⎞⎟ b) x = ⎝ ⎠ = 11,64hs b V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 12 Exercício proposto 3 V – Método dos Mínimos Quadrados – Cálculo Numérico – Prof. Dr. Sergio Pilling 13