Não foi possível enviar o arquivo. Será algum problema com as permissões?
Diferenças
Aqui você vê as diferenças entre duas revisões dessa página.
Próxima revisão | Revisão anterior | ||
disciplinas:ce709-2010:atividades [2010/09/28 16:34] paulojus criada |
disciplinas: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> | ||
+ | |||