Diferenças

Aqui você vê as diferenças entre duas revisões dessa página.

Link para esta página de comparações

Próxima revisão
Revisão anterior
disciplinas:ce709-2010:atividades [2010/09/28 16:34] – criada paulojusdisciplinas:ce709-2010:atividades [2010/10/20 20:01] (atual) walmes
Linha 1: Linha 1:
 ====== Atividades e Materiais ====== ====== Atividades e Materiais ======
 +
 +=====Modelos não lineares=====
 +
 +Saudações. Está previsto para a disciplina de estudos dirigidos em estatística a abordagem do tópico Modelos Não Lineares. Para início de discussão os seguintes materiais devem ser lidos:
 +
 +  * Capítulo //8 Non-Linear and Smooth Regression// do //Modern Applied Statistics with S-plus// ([[http://www.leg.ufpr.br/~walmes/docs/Modern%20Applied%20Statistical%20Methods%20with%20S%20plus.pdf |MASS]]);
 +  * Apêndice do John Fox sobre regressão não linear com o R ([[http://cran.r-project.org/doc/contrib/Fox-Companion/appendix-nonlinear-regression.pdf |Nonlinear Regression and Nonlinear Least Squares]]);
 +  * Exemplo de ajuste de modelo não linear com o R da página do Prof Paulo Justiniano ([[http://leg.ufpr.br/~paulojus/embrapa/Rembrapa/Rembrapase31.html#x33-19000031 |Ajuste de modelos não lineares]]);
 +  * Fazer o estudo dos dados abaixo. Propor um modelo não linear para os dados. Fazer o ajuste no R com as informações e rotinas vistas no materiais acima.
 +
 +<code R>
 +#-----------------------------------------------------------
 +# potássio liberado acumulado para o esterco de codorna
 +
 +klib <- c(51.03, 57.76, 26.60, 60.65, 87.07, 64.67,
 +          91.28, 105.22, 72.74, 81.88, 97.62, 90.14,
 +          89.88, 113.22, 90.91, 115.39, 112.63, 87.51,
 +          104.69, 120.58, 114.32, 130.07, 117.65, 111.69,
 +          128.54, 126.88, 127.00, 134.17, 149.66, 118.25,
 +          132.67, 154.48, 129.11, 151.83, 147.66, 127.30)
 +
 +#-----------------------------------------------------------
 +# tempo em que foram feitas as coletas
 +
 +tempo <- rep(c(15, 30, 45, 60, 75, 90,
 +               120, 150, 180, 210, 240, 270), each=3)
 +liber <- data.frame(tempo, klib)
 +
 +require(lattice)
 +
 +#-----------------------------------------------------------
 +# previa gráfica dos dados
 +
 +xyplot(klib~tempo, data=liber,
 +       type=c("p", "smooth"), col=1,
 +       xlab="Período de incubação (dias)",
 +       ylab="Potássio liberado acumulado (mg/kg de solo)")
 +
 +#-----------------------------------------------------------
 +# conjunto de dados com as médias das repetições e prévia
 +
 +lmedio <- data.frame(tempo=unique(liber$tempo),
 +                     kmedio=tapply(liber$klib,
 +                       liber$tempo, mean))
 +
 +xyplot(kmedio~tempo, data=lmedio,
 +       type=c("p", "smooth"), col=1,
 +       xlab="Período de incubação (dias)",
 +       ylab="Potássio liberado acumulado (mg/kg de solo)")
 +
 +#-----------------------------------------------------------
 +# daqui em diante é com você...
 +
 +#-----------------------------------------------------------
 +# exemplo de uso da optim()
 +
 +x <- 1:9
 +A <- 5
 +B <- 1
 +y <- A*x/(B+x)+rnorm(x,0,0.1)
 +plot(y~x)
 +curve(A*x/(B+x), add=TRUE)
 +
 +#-----------------------------------------------------------
 +# definição da função objetivo
 +
 +fun.objetivo <- function(theta, y, x){
 +  sum((y-theta[1]*x/(theta[2]+x))^2)
 +}
 +
 +#-----------------------------------------------------------
 +# escolha de valores iniciais
 +
 +start <- c(3,0.5)
 +
 +#-----------------------------------------------------------
 +# optimização da função objetivo
 +
 +opt <- optim(start, fun.objetivo, y=y, x=x)
 +opt
 +
 +curve(opt$par[1]*x/(opt$par[2]+x), add=TRUE, col=2)
 +
 +#-----------------------------------------------------------
 +# usando outra função objetivo
 +
 +fun.objetivo <- function(theta, y, x){
 +  n <- length(y)
 +  -(-n/2*log(2*pi)-n/2*log(theta[3])-
 +    sum((y-theta[1]*x/(theta[2]+x))^2/(2*theta[3])))
 +}
 +
 +#-----------------------------------------------------------
 +# os chutes
 +
 +start <- c(3,0.5,0.1)
 +
 +#-----------------------------------------------------------
 +# optimização
 +
 +opt <- optim(start, fun.objetivo, y=y, x=x)
 +opt
 +
 +curve(opt$par[1]*x/(opt$par[2]+x), add=TRUE, col=3)
 +#-----------------------------------------------------------
 +
 +#------------------------------------------------------------------------------------------
 +
 +library(gWidgetsRGtk2)
 +
 +da <- data.frame(x=1:20)
 +da$y <- 10*da$x/(3+da$x)+rnorm(da$x,0,0.2)
 +plot(y~x, data=da)
 +
 +#------------------------------------------------------------------------------------------
 +
 +limits <- list(A=c(0,20), B=c(0,6), n=c(0,2))
 +
 +plotMM <- function(...){
 +  plot(y~x, data=da)
 +  curve(svalue(A)*x^svalue(n)/(svalue(B)+x), add=TRUE)
 +}
 +
 +  w <- gwindow("Slider and spinbox example")
 +
 +  tbl = glayout(cont=w)
 +  for(i in 1:length(limits)){
 +    tbl[i,1] <- paste("Slide to adjuste parameter", names(limits)[i])
 +    tbl[i,2, expand=TRUE] <- (assign(names(limits)[i],
 +               gslider(from=limits[[i]][1],
 +                       to=limits[[i]][2],
 +                       by=diff(limits[[i]])/20,
 +                       value=mean(limits[[i]]),
 +                       container=tbl, handler=plotMM)))
 +  }
 +
 +plotMM()
 +
 +#------------------------------------------------------------------------------------------
 +</code>
 +
  

QR Code
QR Code disciplinas:ce709-2010:atividades (generated for current page)