[R-br] MÉTODO DE GAUSS-NEWTON
Tiago Souza Marçal
tiagosouzamarcal em hotmail.com
Sábado Janeiro 19 18:00:47 BRST 2013
Boa tarde colegas,
Estou tentando reproduzir o algoritmo de GAUSS-NEWTOM em linguagem r, entretanto não estou conseguindo automatizar os comandos.
Na parte da automatização eu estou seguindo o script que o Walmes postou no site <http://ridiculas.wordpress.com/?s=newton> sobre o metodo de NEWTON-RAPHSON.
A titilo de exemplo estou utilizando um ajuste para o modelo y=a*exp(b*x) conforme segue a baixo:
CMR:
x<-c(2,5,7,10,14,19,26,31,34,38,45,52,53,60,65)
y<-c(54,50,45,37,35,25,20,16,18,13,8,11,8,4,6)
a<-56.6646 #parâmetro obtido obtido pela linearização da função y=a*exp(b*x)
b<-(-0.03797) #parâmetro obtido obtido pela linearização da função y=a*exp(b*x)
# a função fx fornece as estimativas de (y~x) dados os coeficientes a e b.
fx<-function(xi,yi,a,b,i,n){resp<-numeric(0)for(i in i:n){t<-(a*exp(b*xi[i]))resp<-c(resp,t) }resp}
# a função dfa fornece as estimativas da primeira derivada (dfx/da ~ x) dado o coeficientes b.
dfa<-function(xi,yi,b,i,n){resp<-numeric(0)for(i in i:n){ t<-exp(b*xi[i]) resp<-c(resp,t) } resp }
# a função dfb fornece as estimativas da primeira derivada (dfx/db~x) dados os coeficientes a e b.
dfb<-function(xi,yi,a,b,i,n){resp1<-numeric(0)for(i in i:n){ t1<-a*xi[i]*exp(b*xi[i]) resp1<-c(resp1,t1) } resp1 }
# A matrix betas contem os chutes para a e b
betas<-matrix(c(a,b),nc=1)
# bcor fornece os parâmetros a e b corrigidos para a primeira interação. Estes parâmetros serão substituídos no lugar dos valores de a e b iniciais e o processo sera repetido.
bcor<-(betas)+((solve((t(matrix(c((dfa(dia,diag,b,1,15)),(dfb(dia,diag,a,b,1,15))),nc=2))%*%(matrix(c((dfa(dia,diag,b,1,15)),(dfb(dia,diag,a,b,1,15))),nc=2)))))%*%((t(matrix(c((dfa(dia,diag,b,1,15)),(dfb(dia,diag,a,b,1,15))),nc=2))%*%(matrix(c(diag-(fx(dia,diag,a,b,1,15))),nc=1)))))
# a função Sqerro fornece a soma de quadrados dos desvios
Sqerro<-sum((diag-fx(dia,diag,a,b,1,15))^2)
O problema é que eu não estou conseguindo usar os parâmetros corrigidos para fazer a nova interação:
Ex:
bcor[i+1]<-bcor[i]+betas[i] # este comando aproximaria mais vezes os parâmetros.
dif(Sqerro[i+1]-Sqerro[i]) #Este seria o critério de parada quando a variação no Sqerro fosse insignificante (10^(-6) por exemplo)
Agradeço desde já pela ajuda e peço desculpas se não fui muito claro.
Att.
Tiago.
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20130119/9405ecd9/attachment.html>
Mais detalhes sobre a lista de discussão R-br