Resumos
O ajustamento polinomial a dados experimentais pelo método de quadrados mínimos é muito usado, e aqui estudamos a introdução de uma simples modificação, que, em certos casos, é mais adequada, isto é, condicionamos o polinômio a conter um dado ponto. Incluímos uma implementação do correspondente algoritmo em MATLAB.
ajustamento; quadrados mínimos; condicionamento
Since the least-squares polynomial fit is used very often, we study here a simple modification that is more suitable for some cases, that is, we constrain the polynomial to pass through a given point, besides satisfying the least-squares property. We include an implementation of the corresponding algorithm in MATLAB.
data fitting; least squares; constraint
NOTAS E DISCUSSÃO
Polinômio dos quadrados mínimos condicionado
The constrained least squares polynomial
Oclide J. Dotto1 1 E-mail: ojdotto@ucs.br. ; Adalberto A. Dornelles Filho
Departamento de Matemática e Estatística, Universidade de Caxias do Sul, Caxias do Sul, RS, Brasil
RESUMO
O ajustamento polinomial a dados experimentais pelo método de quadrados mínimos é muito usado, e aqui estudamos a introdução de uma simples modificação, que, em certos casos, é mais adequada, isto é, condicionamos o polinômio a conter um dado ponto. Incluímos uma implementação do correspondente algoritmo em MATLAB.
Palavras-chave: ajustamento, quadrados mínimos, condicionamento.
ABSTRACT
Since the least-squares polynomial fit is used very often, we study here a simple modification that is more suitable for some cases, that is, we constrain the polynomial to pass through a given point, besides satisfying the least-squares property. We include an implementation of the corresponding algorithm in MATLAB.
Keywords: data fitting, least squares, constraint.
1. Introdução
É popular o ajustamento polinomial discreto dos quadrados mínimos (q. m.). Para evocar o processo desse ajustamento, consideremos os dados da Tabela 1, onde a primeira linha forma uma seqüência finita crescente.
Suponhamos que queiramos ajustar a esses dados o polinômio
de ordem k < n dos q. m. A melhor maneira de obter a solução é o uso da álgebra linear como segue. Começamos gerando a matriz de planejamento (design matrix [1])
de ordem (n+1)×(k+1), determinada pelo vetor x = [ x0x1 ... xn ]t. A solução única (porque X é não-singular) da equação normal
é o vetor a dos coeficientes do polinômio dos q. m. (1).
Do ponto de vista dos métodos de resolução numérica de (3), o uso da matriz de planejamento pode levar a imprecisões uma vez que X se torna mal-condicionada quando sua ordem é alta. No entanto, isso não é relevante na maioria dos casos práticos, onde dificilmente a ordem do polinômio ajustado é maior que 3.
O desenvolvimento mais completo da teoria e solução numérica dos quadrados mínimos foi feita por Åke Björck em [2], mas pode ser encontrado em muitos bons livros de álgebra linear, como [3, 4, 5]. Uma introdução rápida sobre a solução dos q. m. encontra-se em [6].
Façamos aqui uma pequena digressão útil. Na estimativa da solução Xa = y mediante o método dos q. m., é conveniente às vezes ponderar essa estimativa, isto é, usar e equação normal ponderada
em vez da equação normal ordinária (3), escolhendo, para matriz de ponderação W, uma conveniente matriz simétrica definida positiva, que minimize o vetor residual r : = Xa - y no sentido dos q. m. associados ao produto interno definido por W, o que significa simplesmente minimizar rtWr. É conhecido [7] que a escolha ótima é a inversa da matriz V = [ vij ] das variâncias-covariâncias dos vetores aleatórios das observações, yi = [ y1iy2i ... ymi]t, i = 0, 1, ..., n, isto é, W = V-1 (cf. também [2, 8]). Ainda, caso os desvios possam ser supostos não-correlacionados, todas as covariâncias são nulas e V se reduz a uma matriz diagonal, cujos elementos diagonais são as variâncias, que medem a dispersão dos desvios em torno da média, e, nesse caso, para inverter V, basta inverter esses elementos diagonais (que são não-nulos). O efeito do uso de V-1 para ponderação é terem menor peso os desvios maiores.
Voltemos agora ao assunto central do artigo. A curva que representa geometricamente o polinômio dos q. m. não passa pelos pontos definidos pelos dados na Tabela 1, a não ser que coincida com a do polinômio de interpolação, o que ocorre quando k = n. No entanto, em certos casos, pela natureza do problema, sabemos de antemão que deve passar por determinado ponto, mesmo que esse ponto não figure entre os dados.
A um polinômio dos q. m., forçado a passar por um ponto P0 escolhido previamente, chamamos de polinômio dos q. m. condicionado, ou polinômio dos q. m. ancorado no ponto P0.
2. Reta dos quadrados mínimos ancorada na origem
A lei de Hooke estabelece que, quando uma força F é aplicada a uma mola homogênea, o alongamento x da mola varia linearmente com essa força, dentro de limites, isto é,
onde a é dita a constante de Hooke da mola. Para determinar a, executamos experimentos que fornecem uma tabela de duas linhas com os valores de F dados e os valores de x resultantes. Em seguida ajustamos a reta dos q. m. aos dados e a declividade dessa reta é uma estimativa da constante de Hooke da mola. Essa reta quase certamente não passa pela origem, mesmo que a origem se inclua nos dados. No entanto, idealmente, deveria passar pela origem, pois o alongamento é nulo para uma força aplicada nula.
Para forçar a reta dos q. m. a passar pela origem, basta achar a solução dos q. m. do sistema linear
Uma relação de proporcionalidade direta entre duas grandezas, como em (5), ocorre com freqüência. Outros exemplos são: a intensidade i e a tensão v de uma corrente elétrica aplicada a um resistor, a força F e aceleração resultante a aplicada a um móvel, etc.
Exemplo 1
Vejamos um caso específico. A Tabela 2 mostra os valores do alongamento x (em centímetros) de uma mola vs. força de tração F (em Newtons) aplicada, ambos em módulo.
Aqui a reta dos q. m. tem por equação
e a reta dos q. m. ancorada na origem
A Fig. 1 mostra ambas as retas. A primeira não contém a origem. Vemos que as declividades, em unidades coerentes, não coincidem, embora sejam próximas, e é razoável dizer que o valor a = 0,76 é uma estimativa melhor que a = 0,78 para a constante da mola.
Para a reta dos q. m. ancorada na origem, a equação matricial (3) pode ser reescrita como
De (9) deduzimos que a declividade a dessa reta é dada por
Em laboratórios de ensino, em vez da fórmula (10), é usada comumente a fórmula
que expressa a média aritmética a das declividades yi/xi das retas, cada uma das quais contendo a origem e um ponto (xi,yi). A expressão (11) é bastante intuitiva e produz um resultado levemente diferente (a = 0,73, no caso) e, por outro lado, parece não ter nenhuma relação com (10). De fato, elas não são equivalentes. Mas é interessante observar o que segue.
Primeiro, adotando a equação normal ponderada (4) e a inversa V-1 da matriz das variâncias-covariâncias para matriz de ponderação, a fórmula (9) passa a escrever-se
Segundo, suponhamos que, para cada i Î {0, 1, 2, ...,n }, tenhamos feito um conjunto de m observações y1i, y2i, ..., ymi; então podemos falar da matriz das variâncias-covariâncias dos vetores aleatórios dessas observações, yi = [ y1i y2i ... ymi]t, i = 0, 1, 2, ..., n. Usando a fórmula normal (12), ponderada pela inversa da matriz V = [ vij ] das variâncias-covariâncias, ou ponderando as declividades yi/xi em (11) pelos inversos das correspondentes variâncias dessas declividades (ao invés de tomar a média aritmética), obtemos por ambos os caminhos a fórmula
onde pusemos vi : = vii (variâncias) e admitimos vij = 0, se i ¹ j (covariâncias nulas). Com relação a ponderação em (13), lembramos que, por uma propriedade da variância, temos:
var ,
e o denominador principal em (13) é a soma dos inversos dessas variâncias.
Notamos que a fórmula (10) é recuperada de (13) supondo v0 = v1 = ... = vn, e, com uma simplificação maior, que a fórmula (11) segue de (13), usando as correspondentes variâncias do vetor das declividades
todas iguais a 1/(n+1), isto é, para i = 0, 1, 2 ..., n. Isso explica porque as fórmulas (10) e (11) não são equivalentes.
3. Polinômio dos quadrados mínimos condicionado
Em vez de ancorar na origem a reta dos q. m., podemos ancorá-la num ponto qualquer. Mas tratemos do caso geral: mostremos como ancorar o polinômio de ordem k dos q. m. num ponto escolhido (x0, y0). O procedimento é feito em três etapas:
1. operamos uma translação de eixos de maneira que a nova origem seja o ponto (x0, y0); para isso basta subtrair x0 das abscissas e y0 das ordenadas dos dados;
2. ajustamos aos dados transformados o polinômio q dos q. m. ancorado na origem, como no Exemplo 1;
3. fazemos, por último, a transformação inversa da executada na etapa 1, obtendo o polinômio final dos q. m., p(x) = q(x - x0), ancorado em (x0, y0).
Suponhamos que os dados na Tabela 1 já tenham sofrido a alteração descrita na etapa 1. Para operacionalizar a etapa 2, suprimimos a última coluna da matriz de planejamento (2), obtendo uma matriz que indicaremos com U, e achamos a solução a dos q. m. do sistema linear
O vetor-solução a será formado pelos coeficientes do polinômio q ancorado na origem, descrito na etapa 2. Para o caso de uma reta, a equação matricial (15) torna-se uma equação escalar, da forma (9), onde x e y são os vetores dos dados.
No apêndiceapêndice implementamos esse algoritmo em MATLAB, com o nome polqmcond. O programa produz também a correspondente visualização geométrica.
Exemplo 2
Consideremos a integral elíptica completa
Mediante integração numérica, obtemos os valores de E(x) (com x dado em graus) e construímos a Tabela 3.
Se olharmos para a distribuição desses dados, veremos que um ajustamento parabólico é adequado. Além disso, é óbvio que E(0) = p/2. Então é indicado ancorar a parábola dos q. m. no ponto (0, p/2). Utilizando o programa polqmcond, obtemos essa parábola
A parábola ''livre'' dos q. m. é dada por
A Fig. 2 exibe ambas as parábolas. Visualmente elas se superpõem em quase todo o trajeto.
Para compararmos os dois ajustamentos, os testamos nos pontos x = 2, 12, 17, 27. Obtivemos a Tabela 4 na qual a segunda linha contém os valores de E(x) obtidos por integração numérica aqui considerados exatos, a terceira linha os valores de P(x) do ajustamento dos q. m. ancorado, e a quarta linha os valores de Q(x) do ajustamento livre.
Vemos que, nesse exemplo, o ajustamento ancorado é melhor para os pontos próximos ao ponto ancorador, como era esperado, e longe dele os ajustamentos tendem a ser equivalentes. Claramente, para um ajustamento linear, por exemplo, a situação longe da âncora é outra, como mostra a Fig. 1.
4. Comentários finais
Nos laboratórios de ensino, é bastante comum realizar experimentos para obter o coeficiente de proporcionalidade entre duas grandezas diretamente proporcionais. Se for usado o recurso do ajuste linear dos q. m., disponível em algumas calculadoras científicas, dificilmente obteremos uma reta que passe exatamente pela origem. Uma elegante solução para essa dificuldade é usar a fórmula (10) que decorre do ajustamento dos q. m. ancorado na origem.
Mais geral, mostramos que é possível e relativamente fácil, ao ajustarmos um polinômio dos quadrados mínimos, condicioná-lo a conter um ponto de interesse. Poder-se-ia tentar desenvolver um algoritmo que o condicione a conter mais de um ponto, contanto que ainda retenha um número de graus de liberdade suficientemente elevado. Por exemplo, um caso extremo oposto é condicionar uma reta dos q. m. a conter dois pontos, o que corresponde a ela ficar com zero graus de liberdade e equivale ao problema trivial de achar a equação da reta que passa por esses dois pontos. Não temos certeza se há interesse prático nessa generalização de ancoragem. De qualquer maneira, parece-nos que o caso interessante, exeqüível e útil, é a ancoragem em um ponto apenas.
Agradecimento
Os autores agradecem ao parecerista a gentil apreciação do artigo, assim como as sugestões.
Recebido em 20/2/2006; Revisado em 15/5/2006; Aceito em 29/5/2006
- [1] Richard A. Johnson and Dean W. Wichern, Applied Multivariate Analysis (Prentice Hall, Nova York, 2002)
- [2] Ĺke Björck, Numerical Methods for Least Squares Problems (SIAM, Philadelphia, 1966).
- [3] James W. Demmel, Applied Numerical Linear Algebra (SIAM, Philadelphia, 1997).
- [4] Gene H. Golub and Charles F. Van Loan, Matrix Computations (The John Hopkins University Press, London, 1996).
- [5] Lloyd N. Trefethen and David Bau III, Numerical Linear Algebra (SIAM, Philadelphia, 1997).
- [6] Mário Barone Júnior, Álgebra Linear (IME-USP, Săo Paulo, 2005), 3Ş ed.
- [7] Gilbert Strang and Kay Borre, Linear Algebra, Geodesy and GPS (Wellesley-Cambridge Press, Wellesley, 1997).
- [8] José Henrique Vuolo, Fundamentos da Teoria de Erros (Ed. Edgard Blücher Ltda., Săo Paulo, 2005), 2Ş ed.
apêndice
![](https://minio.scielo.br/documentstore/1806-9126/kdJTs4xk85B7s4QCLL8Fwwy/672d2d0524c17706f5860950dc9dd31543977e44.gif)
Datas de Publicação
-
Publicação nesta coleção
06 Nov 2006 -
Data do Fascículo
2006
Histórico
-
Aceito
29 Maio 2006 -
Revisado
15 Maio 2006 -
Recebido
20 Fev 2006