Este tipo de problema é muito frequente em ciências experimentais. Para ter um exemplo em mente durante a discussão que será feita mais adiante, considere as seguintes informações:
Segundo dados disponibilizados pelo IBGE, o Produto Interno Bruto per capita na cidade de São Paulo, no período de 2002 a 2005, foi (em reais): 17734, 19669, 20943 e 24083, respectivamente.
Com base nessas informações, como poderia ser feita uma previsão do valor correspondente ao ano seguinte (2006)?
Uma escolha possível seria supor que a cada ano o PIB aumenta uma aproximadamente constante, ou seja, usar um modelo linear para obter tal estimativa (que possivelmente será bem grosseira). Intuitivamente, bastaria analisar os dados disponíveis e a partir deles deduzir qual é o aumento que ocorre a cada ano. Depois, a previsão para 2006 seria aproximadamente igual à de 2005 somada com aquele aumento anual.
Esta idéia poderia até funcionar para o caso deste exemplo, mas o que fazer se a quantidade de dados disponíveis sobre algum fenômeno (ou alguma situação) for significativamente maior?
A melhor escolha, sem dúvida, é fazer uso de um computador para obter o modelo que melhor descreve o "comportamento" dos dados experimentais.
Em geral, os problemas de mínimos quadrados consistem em identificar os valores de determinados parâmetros, de modo que se satisfaçam certas equações
. No contexto do exemplo anterior, se procura um modelo linear para os dados, ou seja, uma função
que os descreva da melhor forma possível. Assim, os parâmetros a considerar são
e
. Os valores ideais para essas variáveis seriam aqueles que verificassem as seguintes equações:
![{\displaystyle \left\{{\begin{matrix}2002m+n&=&17734\\2003m+n&=&19669\\2004m+n&=&20943\\2005m+n&=&24083\end{matrix}}\right.}](https://wikimedia.org/api/rest_v1/media/math/render/svg/89b83cf684f165606f13d66aa29b1043b01d39df)
Note que a partir dessas equações poderiam ser definidas as funções
como:
![{\displaystyle h_{1}(m,n)=2002m+n-17734}](https://wikimedia.org/api/rest_v1/media/math/render/svg/217574eb8c318b7534cde301b0299d18c4c5408f)
![{\displaystyle h_{2}(m,n)=2003m+n-19669}](https://wikimedia.org/api/rest_v1/media/math/render/svg/6e2ff59dfb804d5ff4bb332abb86c33fd17e3d58)
![{\displaystyle h_{3}(m,n)=2004m+n-20943}](https://wikimedia.org/api/rest_v1/media/math/render/svg/070d47487473e2db7434b36cd604a186b4545e4a)
![{\displaystyle h_{4}(m,n)=2005m+n-24083}](https://wikimedia.org/api/rest_v1/media/math/render/svg/dba5fff32ecc3bb5c99333acc521319b20a651db)
É de se esperar que o sistema de equações obtido a pouco não admitirá uma solução exata, pois se tem mais equações do que variáveis.
Isso geralmente acontece, pois é comum haver uma quantidade
de equações bem maior que o número
de parâmetros a identificar. Em particular, quando todas as equações
são lineares em suas
variáveis, dificilmente existirá uma solução exata para o sistema linear resultante, pois este terá mais equações do que incógnitas (como no exemplo). Em geral, não é possível encontrar parâmetros que satisfaçam exatamente todas as equações. Por isso, costuma-se tentar identificar os parâmetros que "melhor se aproximam" de uma solução exata, em algum sentido.
Uma forma de obter uma solução aproximada (uma "quase-solução") resulta da seguinte observação: o valor de cada função
em uma solução exata
deveria ser zero. Se tal exigência é restritiva demais, e com ela não é possível encontrar qualquer solução, uma possibilidade seria exigir um pouco menos. Por exemplo, poderia ser exigido apenas que o valor de
, para
seja, em geral, pequeno. Uma das formas de capturar essa idéia em termos mais precisos é dizer que se pretende minimizar a soma dos quadrados dos valores de cada
. Em símbolos, o problema passaria a ser:
![{\displaystyle \min _{x\in \mathbb {R} ^{n}}\sum _{i=1}^{p}\left[h_{i}(x)\right]^{2}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/31c65c0dccda26b6c814e23c07fe3feb3edbe0eb)
Neste caso, para cada índice
, a função
é afim linear, ou seja:
![{\displaystyle h_{i}:\mathbb {R} ^{n}\mapsto \mathbb {R} }](https://wikimedia.org/api/rest_v1/media/math/render/svg/2c2ae945c20024e5943cde1db482291f747c2ff4)
![{\displaystyle h_{i}(x)=a_{i}^{\top }x+b_{i}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1bb619bb5cb1d04cc8f1b018e44e8e954f410bf7)
onde
e
para cada
. Deste modo, pode-se definir uma função
como
de modo que
![{\displaystyle H(x)={\begin{bmatrix}h_{1}(x)\\\vdots \\h_{p}(x)\end{bmatrix}}={\begin{bmatrix}a_{1}^{\top }x+b_{1}\\\vdots \\a_{p}^{\top }x+b_{p}\end{bmatrix}}={\begin{bmatrix}a_{1}^{\top }\\\vdots \\a_{p}^{\top }\end{bmatrix}}x+{\begin{bmatrix}b_{1}\\\vdots \\b_{p}\end{bmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/37c847f4d35cfca0b84954dab535883ad7099744)
Motivando a introdução da seguinte notação:
![{\displaystyle A={\begin{bmatrix}a_{1}^{\top }\\\vdots \\a_{p}^{\top }\end{bmatrix}},\quad b={\begin{bmatrix}b_{1}\\\vdots \\b_{p}\end{bmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a1739fe09633c39da0e9ce2ddc5e7ada87f75d86)
Assim,
![{\displaystyle H(x)=Ax+b}](https://wikimedia.org/api/rest_v1/media/math/render/svg/75483595ea63aec31469b1faa852262b536c1292)
Logo, buscar uma solução exata é o mesmo que procurar uma solução para o seguinte sistema linear:
![{\displaystyle Ax=-b}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ae722182aaad594757eab2f93499efe4d5103ca3)
E uma solução aproximada poderia ser buscada a partir do seguinte problema de minimização:
![{\displaystyle \min _{x\in \mathbb {R} ^{n}}{\frac {1}{2}}\|Ax+b\|^{2}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c8fc5656ef607ae444e69110b9c03c2031de660f)
- Exercício
Verifique que se houver uma solução exata, ela também será um minimizador de
. A recíproca é verdadeira?
Analisando a função objetivo do problema de minimização anterior, tem-se:
![{\displaystyle f(x)={\frac {1}{2}}\|Ax+b\|^{2}={\frac {1}{2}}\langle Ax+b,Ax+b\rangle ={\frac {1}{2}}x^{\top }A^{\top }Ax+b^{\top }Ax+{\frac {1}{2}}b^{\top }b}](https://wikimedia.org/api/rest_v1/media/math/render/svg/214027ad551e4a1fdc1245adeacd01e428dbcc09)
Logo, como
é simétrica e semi-definida positiva, tem-se
convexa. Isso implica que a condição necessária de primeira ordem é também suficiente. Assim, qualquer ponto
tal que
é solução do problema aproximado.
Calculando o gradiente da função objetivo tem-se:
![{\displaystyle 0=\nabla f(x)=A^{\top }Ax+A^{\top }b}](https://wikimedia.org/api/rest_v1/media/math/render/svg/4723c4e0d4f55f1473659a74da11bf1c1634ab9f)
Deste modo, a solução do problema é obtida resolvendo o sistema
![{\displaystyle A^{\top }Ax=-A^{\top }b}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a6ce0a24827cf127e544e387170165db5d615bcf)
Observe também que isso implica
, ou seja,
.
Considere dado um conjunto de pontos do plano
, por exemplo
, representando dados obtidos experimentalmente.
Perguntas:
- 1. Qual é a função afim linear que melhor se aproxima dos dados experimentais?
- 2. Qual é a função quadrática que melhor se aproxima dos dados experimentais?
Resolução
|
As funções quadráticas possuem a forma . Então os parâmetros a identificar são , e .
Pode-se tomar .
Procedendo como antes, tem-se
![{\displaystyle H(x)=Ax+c}](https://wikimedia.org/api/rest_v1/media/math/render/svg/3143b50a61dca788d31afe9a0b4765ecad2327fc)
onde
![{\displaystyle A={\begin{bmatrix}x_{1}^{2}&x_{1}&1\\\vdots &\vdots &\vdots \\x_{p}^{2}&x_{p}&1\end{bmatrix}},\quad x={\begin{bmatrix}a\\b\\c\end{bmatrix}},\quad c={\begin{bmatrix}-y_{1}\\\vdots \\-y_{p}\end{bmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/4add2157967b6c4e4e3169192d7d030e53926bd8)
|
- Exercício
Explorar o exemplo anterior com dados concretos, implementando um dos três métodos de gradientes conjugados.
Resolução
|
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a neste módulo.
|
- Observações
Conforme se aumenta o grau do polinômio que faz a aproximação dos dados, as colunas de
têm elementos elevados a potências cada vez maiores, fazendo com que os autovalores de
sejam cada vez mais dispersos. Com isso,
torna-se mal condicionada.
- Exercício
É possível, usando o problema de mínimos quadrados, verificar o postulado de Euclides que diz que "por dois pontos distintos passa uma única reta"?
Resolução
|
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a neste módulo.
|
- Exercício
Usando o problema de mínimos quadrados, é possível inferir que "por três pontos distintos passa uma única curva quadrática"?
Resolução
|
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a neste módulo.
|
Para esse tipo de problemas, há dois métodos:
- Gauss-Newton
- Levemberg-Marquardt
Ambos são métodos do tipo Newton. Então, para entender cada um deles é preciso entender o Método de Newton.
Para entender a essência do método de Newton, primeiro considere que o problema a ser resolvido é
![{\displaystyle (P)\left\{{\begin{matrix}\min f(x)\\x\in \mathbb {R} ^{n}\end{matrix}}\right.}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e014df26517fc35eb8ef61b57aa756f0f29d4d10)
sendo
, e portanto
de classe
. A idéia de Newton é usar o desenvolvimento até segunda ordem da série de Taylor da função
em cada ponto iterado. Isto é, se o iterado é
, então:
![{\displaystyle Q(x)=f({\bar {x}})+\nabla f({\bar {x}})^{\top }(x-{\bar {x}})+{\frac {1}{2}}(x-{\bar {x}})\nabla ^{2}f({\bar {x}})(x-{\bar {x}})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8cff760d39daede80e47442b2fe5c204b1282076)
Então a condição de Newton é que em cada iteração a Hessiana
deve ser definida positiva.
Calculando o gradiente da função
, segue:
![{\displaystyle \nabla Q(x)=\nabla f({\bar {x}})+\nabla ^{2}f({\bar {x}})(x-{\bar {x}})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d1e59878c8e33aa205949b64ac04f853013192ea)
Se
é o (único) minimizador de
, então
![{\displaystyle 0=\nabla Q(x^{*})=\nabla f({\bar {x}})+\nabla ^{2}f({\bar {x}})(x^{*}-{\bar {x}})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/3bef7c25aa7b71b2838bee9f8e9961909e10c44c)
donde
![{\displaystyle \nabla ^{2}f({\bar {x}})(x^{*}-{\bar {x}})=-\nabla f({\bar {x}})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/af095a189f1dad5ca0a20764945b737665d6a440)
Sendo
definida positiva, tal matriz é também inversível. Portanto:
![{\displaystyle x^{*}={\bar {x}}-\left[\nabla ^{2}f({\bar {x}})\right]^{-1}\nabla f({\bar {x}})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/65d5e7ce18c97bc0920de31d09c5dd24d6d1eda5)
Assim, pode-se usar a seguinte iteração:
![{\displaystyle x_{k+1}=x_{k}-\left[\nabla ^{2}f(x_{k})\right]^{-1}\nabla f(x_{k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/70ece61409bcdab05cb09433ee4dfcae51783720)
Início: Tome
Se
, pare:
é ponto crítico.
Senão, Calcule
, solução de
Faça
e
Iteração: Se
, pare:
é ponto crítico.
Senão, calcular
, solução de
Faça
e
Voltando ao problema original, de mínimos quadrados, se tinha:
![{\displaystyle f(x)={\frac {1}{2}}\sum _{i=1}^{p}\left[h_{i}(x)\right]^{2}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/0ebe363eac1b7b7878a343733dd696c9fbfd802c)
Calculando o gradiente desta função, resulta:
![{\displaystyle \nabla f(x)=\sum _{i=1}^{p}h_{i}(x)\nabla h_{i}(x)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/2d7df70705cefd7eeb6147b6a09ff70f4c060b79)
Considera-se
definida por
![{\displaystyle H(x)={\begin{bmatrix}h_{1}(x)\\\vdots \\h_{p}(x)\end{bmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/511b9c926f44da12fd33ae7b1b92f96c30b5761a)
Deste modo, a Jacobiana de
verifica:
![{\displaystyle \nabla f(x)=\left[J_{H}(x)\right]H(x)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1981dd53fef487cb06c9112714d7065d7261e380)
pois o produto de uma matriz por um vetor tem como resultado um vetor que é a combinação linear das colunas da matriz, com coeficientes que são as coordenadas do vetor.
Além disso, tem-se
![{\displaystyle \nabla ^{2}f(x)=\sum _{i=1}^{p}h_{i}(x)\nabla ^{2}h_{i}(x)+\sum _{i=1}^{p}\nabla h_{i}(x)\nabla h_{i}(x)^{\top }}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ced03dc89f61cceca6218a52eb560abc948e3621)
Seja
![{\displaystyle B(x)=\sum _{i=1}^{p}\nabla h_{i}(x)\nabla h_{i}(x)^{\top }=J_{H}(x)J_{H}(x)^{\top }}](https://wikimedia.org/api/rest_v1/media/math/render/svg/6ab5676b15da0d245241895d20496b5e908aa5c9)
Sabe-se que uma matriz
é definida positiva se
para qualquer
. Fazendo
, tem-se:
![{\displaystyle x^{\top }\nabla h_{i}(x)\nabla h_{i}(x)^{\top }x=\left[\nabla h_{i}(x)^{\top }x\right]^{2}\geq 0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/debb996bd3a9b547879dac97c408d3255cb28015)
Para que
seja definida positiva, é necessário que
( deve gerar todo o espaço), neste caso, se diz que
é de posto máximo.
Início: Tome
Se
, pare:
é ponto crítico.
Senão, Calcule
, solução de
, onde
Faça
e
Iteração: Se
, pare:
é ponto crítico.
Senão, calcular
, solução de
Faça
e
A idéia de Levemberg-Marquardt foi perturbar a matriz
, considerando
, para algum
pequeno.
Início: Tome
Se
, pare:
é ponto crítico.
Senão, Calcule
, solução de
, onde
Faça
e
Iteração: Se
, pare:
é ponto crítico.
Senão, calcular
, solução de
Faça
e
- Exercício
Enumere as diferenças que existem nos seguintes algoritmos:
- Newton
- Gauss-Newton
- Levenberg-Marquardt
Resolução
|
Newton
Desenvolvendo em série de Taylor e exigindo que a matriz hessiana de seja definida positiva o método de Newton usa .
Gauss-Newton
Em alguns casos a matriz hessiana de pode ser bem complicada. Nesses casos então a idéia do método de Gauss-Newton é usar uma boa aproximação para essa matriz.
Podemos escrever
.
Usando Taylor em temos:
onde é a jacobiana de .
Substituindo na temos
Calculando o gradiente e depois a segunda derivada temos:
e
Denotando , então o método de Gauss-Newton escolhe tal que
.
Levenberg-Marquardt
A idéia do método de Levenberg-Marquardt é pertubar a matriz um pouco, considerando com pequeno.
Portanto o método de Levenberg-Marquardt escolhe tal que
.
|