Mínimos quadrados não lineares - Non-linear least squares

Mínimos quadrados não lineares é a forma de análise de mínimos quadrados usada para ajustar um conjunto de m observações com um modelo que é não linear em n parâmetros desconhecidos ( m  ≥  n ). É usado em algumas formas de regressão não linear . A base do método é aproximar o modelo por um linear e refinar os parâmetros por iterações sucessivas. Existem muitas semelhanças com os mínimos quadrados lineares , mas também algumas diferenças significativas . Em teoria econômica, o método de mínimos quadrados não linear é aplicado em (i) a regressão probit, (ii) regressão de limiar, (iii) regressão suave, (iv) regressão de ligação logística, (v) regressores transformados de Box-Cox ( ) .

Teoria

Considere um conjunto de pontos de dados e uma curva (função de modelo) que além da variável também depende de parâmetros, com É desejado encontrar o vetor de parâmetros de forma que a curva se ajuste melhor aos dados dados no sentido de mínimos quadrados, ou seja, a soma dos quadrados

é minimizado, onde os resíduos (erros de predição na amostra) r i são dados por

para

O valor mínimo de S ocorre quando o gradiente é zero. Como o modelo contém n parâmetros, existem n equações de gradiente:

Em um sistema não linear, as derivadas são funções tanto da variável independente quanto dos parâmetros, portanto, em geral, essas equações de gradiente não têm uma solução fechada. Em vez disso, os valores iniciais devem ser escolhidos para os parâmetros. Em seguida, os parâmetros são refinados iterativamente, ou seja, os valores são obtidos por aproximações sucessivas,

Aqui, k é um número de iteração e o vetor de incrementos é conhecido como vetor de deslocamento. Em cada iteração, o modelo é linearizado por aproximação a uma expansão polinomial de Taylor de primeira ordem sobre

O Jacobiano , J , é uma função das constantes, da variável independente e dos parâmetros, então ele muda de uma iteração para a próxima. Assim, em termos do modelo linearizado,

e os resíduos são dados por

Substituindo essas expressões nas equações de gradiente, elas se tornam

que, no rearranjo, tornam-se n equações lineares simultâneas, as equações normais

As equações normais são escritas em notação matricial como

Essas equações formam a base do algoritmo de Gauss-Newton para um problema de mínimos quadrados não linear.

Observe a convenção de sinais na definição da matriz Jacobiana em termos das derivadas. Fórmulas lineares em podem aparecer com fator de em outros artigos ou na literatura.

Extensão por pesos

Quando as observações não são igualmente confiáveis, uma soma ponderada de quadrados pode ser minimizada,

Cada elemento da matriz de peso diagonal W deve, idealmente, ser igual ao recíproco da variância do erro da medição. As equações normais são, então, mais geralmente,


Interpretação geométrica

Em mínimos quadrados lineares, a função objetivo , S , é uma função quadrática dos parâmetros.

Quando houver apenas um parâmetro, o gráfico de S em relação a esse parâmetro será uma parábola . Com dois ou mais parâmetros, os contornos de S em relação a qualquer par de parâmetros serão elipses concêntricas (assumindo que a matriz de equações normais é definida positiva ). Os valores mínimos dos parâmetros encontram-se no centro das elipses. A geometria da função objetivo geral pode ser descrita como elíptica parabolóide. No NLLSQ a função objetivo é quadrática em relação aos parâmetros apenas em uma região próxima ao seu valor mínimo, onde a série de Taylor truncada é uma boa aproximação do modelo.

Quanto mais os valores dos parâmetros diferem de seus valores ideais, mais os contornos se desviam da forma elíptica. Uma consequência disso é que as estimativas iniciais dos parâmetros devem ser o mais próximo possível de seus (desconhecidos!) Valores ótimos. Ele também explica como a divergência pode ocorrer visto que o algoritmo de Gauss-Newton é convergente apenas quando a função objetivo é aproximadamente quadrática nos parâmetros.

Computação

Estimativas de parâmetros iniciais

Alguns problemas de mau condicionamento e divergência podem ser corrigidos encontrando estimativas de parâmetros iniciais que estão perto dos valores ótimos. Uma boa maneira de fazer isso é por simulação de computador . Os dados observados e calculados são exibidos em uma tela. Os parâmetros do modelo são ajustados manualmente até que a concordância entre os dados observados e calculados seja razoavelmente boa. Embora este seja um julgamento subjetivo, é suficiente encontrar um bom ponto de partida para o refinamento não linear. As estimativas de parâmetros iniciais podem ser criadas usando transformações ou linearizações. Melhor ainda, algoritmos evolutivos, como o algoritmo do funil estocástico, podem levar à bacia convexa de atração que circunda as estimativas dos parâmetros ideais. Algoritmos híbridos que usam randomização e elitismo, seguidos por métodos de Newton têm se mostrado úteis e eficientes computacionalmente.

Solução

Qualquer método entre os descritos abaixo pode ser aplicado para encontrar uma solução.

Critérios de convergência

O critério de senso comum para convergência é que a soma dos quadrados não diminui de uma iteração para a próxima. No entanto, este critério é frequentemente difícil de implementar na prática, por várias razões. Um critério de convergência útil é

O valor 0,0001 é um tanto arbitrário e pode precisar ser alterado. Em particular, pode ser necessário aumentar quando os erros experimentais forem grandes. Um critério alternativo é

Novamente, o valor numérico é um tanto arbitrário; 0,001 é equivalente a especificar que cada parâmetro deve ser refinado para a precisão de 0,1%. Isso é razoável quando é menor do que o maior desvio padrão relativo nos parâmetros.

Cálculo do Jacobiano por aproximação numérica

Existem modelos para os quais é muito difícil ou mesmo impossível derivar expressões analíticas para os elementos do Jacobiano. Então, a aproximação numérica

é obtido pelo cálculo de para e . O incremento,, tamanho deve ser escolhido para que a derivada numérica não esteja sujeita a erro de aproximação por ser muito grande, ou erro de arredondamento por ser muito pequena.

Erros de parâmetro, limites de confiança, resíduos etc.

Algumas informações são fornecidas na seção correspondente na página de mínimos quadrados lineares .

Múltiplos mínimos

Múltiplos mínimos podem ocorrer em uma variedade de circunstâncias, algumas das quais são:

  • Um parâmetro é elevado a uma potência de dois ou mais. Por exemplo, ao ajustar dados a uma curva Lorentziana

onde é a altura, é a posição e é a meia largura a meia altura, existem duas soluções para a meia largura, e que dão o mesmo valor ótimo para a função objetivo.

  • Dois parâmetros podem ser trocados sem alterar o valor do modelo. Um exemplo simples é quando o modelo contém o produto de dois parâmetros, pois dará o mesmo valor que .
  • Um parâmetro está em uma função trigonométrica, como , que tem valores idênticos em . Veja o algoritmo de Levenberg – Marquardt para um exemplo.

Nem todos os mínimos múltiplos têm valores iguais da função objetivo. Mínimos falsos, também conhecidos como mínimos locais, ocorrem quando o valor da função objetivo é maior do que seu valor no chamado mínimo global. Para ter certeza de que o mínimo encontrado é o mínimo global, o refinamento deve ser iniciado com valores iniciais dos parâmetros amplamente diferentes. Quando o mesmo mínimo é encontrado, independentemente do ponto de partida, é provável que seja o mínimo global.

Quando existem vários mínimos, há uma consequência importante: a função objetivo terá um valor máximo em algum lugar entre dois mínimos. A matriz de equações normais não é definida positivamente no máximo na função objetivo, pois o gradiente é zero e não existe uma direção única de descida. O refinamento de um ponto (um conjunto de valores de parâmetro) próximo a um máximo será mal condicionado e deve ser evitado como ponto de partida. Por exemplo, ao ajustar um Lorentziano, a matriz de equações normais não é definida positiva quando a meia largura da banda é zero.

Transformação para um modelo linear

Um modelo não linear às vezes pode ser transformado em linear. Por exemplo, quando o modelo é uma função exponencial simples,

ele pode ser transformado em um modelo linear tomando logaritmos.

Graficamente, isso corresponde a trabalhar em um gráfico semi-log . A soma dos quadrados torna-se

Esse procedimento deve ser evitado, a menos que os erros sejam multiplicativos e com distribuição normal de log, pois pode dar resultados enganosos. Isso vem do fato de que quaisquer que sejam os erros experimentais em y , os erros em log y são diferentes. Portanto, quando a soma transformada dos quadrados é minimizada, resultados diferentes serão obtidos tanto para os valores dos parâmetros quanto para seus desvios-padrão calculados. No entanto, com erros multiplicativos que são distribuídos log-normalmente, este procedimento fornece estimativas de parâmetros imparciais e consistentes.

Outro exemplo é fornecido pela cinética de Michaelis-Menten , usada para determinar dois parâmetros e :

.

O enredo Lineweaver – Burk

de contra é linear nos parâmetros e , mas muito sensível a erros de dados e fortemente inclinado para ajustar os dados em um intervalo particular da variável independente .

Algoritmos

Método de Gauss-Newton

As equações normais

pode ser resolvido pela decomposição de Cholesky , conforme descrito em mínimos quadrados lineares . Os parâmetros são atualizados iterativamente

onde k é um número de iteração. Embora esse método possa ser adequado para modelos simples, ele falhará se ocorrer divergência. Portanto, a proteção contra divergências é essencial.

Corte em turnos

Se ocorrer divergência, um expediente simples é reduzir o comprimento do vetor de deslocamento , por uma fração, f

Por exemplo, o comprimento do vetor de deslocamento pode ser dividido sucessivamente pela metade até que o novo valor da função objetivo seja menor que seu valor na última iteração. A fração, f pode ser otimizada por uma pesquisa de linha . Como cada valor de teste de f requer que a função objetivo seja recalculada, não vale a pena otimizar seu valor com muito rigor.

Ao usar o corte com deslocamento, a direção do vetor de deslocamento permanece inalterada. Isso limita a aplicabilidade do método a situações em que a direção do vetor de deslocamento não é muito diferente do que seria se a função objetivo fosse aproximadamente quadrática nos parâmetros,

Parâmetro Marquardt

Se ocorrer divergência e a direção do vetor de deslocamento estiver tão longe de sua direção "ideal" que o corte de deslocamento não é muito eficaz, ou seja, a fração, f necessária para evitar a divergência é muito pequena, a direção deve ser alterada. Isso pode ser obtido usando o parâmetro Marquardt . Neste método, as equações normais são modificadas

onde é o parâmetro Marquardt e I é uma matriz de identidade. Aumentar o valor de tem o efeito de alterar a direção e o comprimento do vetor de deslocamento. O vetor de deslocamento é girado na direção da descida mais íngreme

quando

é o vetor de descida mais íngreme. Portanto, quando se torna muito grande, o vetor de deslocamento se torna uma pequena fração do vetor de descida mais íngreme.

Várias estratégias foram propostas para a determinação do parâmetro Marquardt. Tal como acontece com o corte de deslocamento, é um desperdício otimizar este parâmetro de forma muito rigorosa. Em vez disso, uma vez que foi encontrado um valor que provoca uma redução no valor da função objetivo, esse valor do parâmetro é transportado para a próxima iteração, reduzido se possível ou aumentado se necessário. Ao reduzir o valor do parâmetro Marquardt, há um valor de corte abaixo do qual é seguro defini-lo como zero, ou seja, continuar com o método Gauss-Newton não modificado. O valor de corte pode ser definido igual ao menor valor singular do Jacobiano. Um limite para este valor é dado por .

Decomposição QR

O mínimo na soma dos quadrados pode ser encontrado por um método que não envolve a formação de equações normais. Os resíduos com o modelo linearizado podem ser escritos como

O Jacobiano está sujeito a uma decomposição ortogonal; a decomposição QR servirá para ilustrar o processo.

onde Q é uma matriz ortogonal e R é uma matriz que é particionada em um bloco,, e um bloco zero. é triangular superior.

O vetor residual é multiplicado à esquerda por .

Isso não tem efeito sobre a soma dos quadrados, pois Q é ortogonal. O valor mínimo de S é obtido quando o bloco superior é zero. Portanto, o vetor de deslocamento é encontrado resolvendo

Essas equações são facilmente resolvidas porque R é triangular superior.

Decomposição de valor singular

Uma variante do método de decomposição ortogonal envolve a decomposição de valor singular , em que R é diagonalizado por outras transformações ortogonais.

onde é ortogonal, é uma matriz diagonal de valores singulares e é a matriz ortogonal dos vetores próprios de ou equivalentemente os vetores singulares corretos de . Neste caso, o vetor de deslocamento é dado por

A relativa simplicidade desta expressão é muito útil na análise teórica de mínimos quadrados não lineares. A aplicação da decomposição de valores singulares é discutida em detalhes em Lawson e Hanson.

Métodos de gradiente

Existem muitos exemplos na literatura científica onde diferentes métodos foram usados ​​para problemas de ajuste de dados não lineares.

A matriz H é conhecida como matriz Hessiana . Embora este modelo tenha melhores propriedades de convergência perto do mínimo, é muito pior quando os parâmetros estão longe de seus valores ótimos. O cálculo do Hessian aumenta a complexidade do algoritmo. Este método não é de uso geral.
  • Método Davidon – Fletcher – Powell . Este método, uma forma de método pseudo-Newton, é semelhante ao anterior, mas calcula o Hessian por aproximações sucessivas, para evitar ter que usar expressões analíticas para as segundas derivadas.
  • Descida mais íngreme . Embora uma redução na soma dos quadrados seja garantida quando o vetor de deslocamento aponta na direção da descida mais íngreme, esse método geralmente tem um desempenho ruim. Quando os valores dos parâmetros estão longe do ótimo, a direção do vetor de descida mais íngreme, que é normal (perpendicular) aos contornos da função objetivo, é muito diferente da direção do vetor de Gauss-Newton. Isso torna a divergência muito mais provável, especialmente porque o mínimo ao longo da direção da descida mais íngreme pode corresponder a uma pequena fração do comprimento do vetor de descida mais íngreme. Quando os contornos da função objetivo são muito excêntricos, devido à alta correlação entre os parâmetros, as iterações de descida mais íngremes, com corte de deslocamento, seguem uma trajetória lenta em zigue-zague em direção ao mínimo.
  • Pesquisa de gradiente conjugado . Este é um método melhorado baseado na descida mais íngreme com boas propriedades de convergência teórica, embora possa falhar em computadores digitais de precisão finita, mesmo quando usado em problemas quadráticos.

Métodos de pesquisa direta

Os métodos de pesquisa direta dependem de avaliações da função objetivo em uma variedade de valores de parâmetros e não usam derivadas de forma alguma. Eles oferecem alternativas para o uso de derivadas numéricas no método de Gauss-Newton e métodos de gradiente.

  • Pesquisa de variável alternada. Cada parâmetro é variado por sua vez, adicionando-se um incremento fixo ou variável a ele e retendo o valor que produz uma redução na soma dos quadrados. O método é simples e eficaz quando os parâmetros não são altamente correlacionados. Possui propriedades de convergência muito pobres, mas pode ser útil para encontrar estimativas de parâmetros iniciais.
  • Pesquisa Nelder – Mead (simplex) . Um simplex neste contexto é um politopo de n  + 1 vértices em n dimensões; um triângulo em um plano, um tetraedro no espaço tridimensional e assim por diante. Cada vértice corresponde a um valor da função objetivo para um determinado conjunto de parâmetros. A forma e o tamanho do simplex são ajustados variando os parâmetros de forma que o valor da função objetivo no vértice mais alto sempre diminua. Embora a soma dos quadrados possa inicialmente diminuir rapidamente, ela pode convergir para um ponto não estacionário em problemas quase-convexos, por um exemplo de MJD Powell.

Descrições mais detalhadas desses e de outros métodos estão disponíveis em Receitas Numéricas , juntamente com código de computador em várias linguagens.

Veja também

Referências

Leitura adicional

  • Kelley, CT (1999). Métodos iterativos para otimização (PDF) . SIAM Frontiers in Applied Mathematics. nº 18. ISBN 0-89871-433-8.
  • Strutz, T. (2016). Ajuste de dados e incerteza: uma introdução prática aos mínimos quadrados ponderados e além (2ª ed.). Springer Vieweg. ISBN 978-3-658-11455-8.