Re: [R-br] Ajustes de modelos não-lineares em gases in vitro

Olá pessoal, eu enviei uma dúvida sobre qual melhor estrutura de correlação residual devo considerar num modelo não linear (logístico bicompartimental) para ajuste de dados com produção de gases in vitro, que possui uma alta dependência residual. Eu organizei melhor o script e os dados, em anexo (novo anexo), em que os dados que estão na planilha, a primeira coluna é o x variável independente representando o tempo. E o y são as mensurações de gases. Então a nova planilha está exatamente as informações no qual preciso obter o ajuste do conjunto. Tentei inserir algumas estruturas de correlação, mas não deu certo. Alguém sabe como devo ajustar esse conjunto de dados (qual melhor estrutura de correlação residual que posso considerar) usando modelos não lineares na função gnls? Adriele #----------------------------------------------------------------------- setwd('C:\\Users\\Adriele\\Desktop') gas<-read.csv("gases_lista_R.csv",sep = ",", header = TRUE) head(gas) attach(gas) #----------------------------------------------------------------------- # Ajuste usando o gnls #Ajuste do modelo exponencial -lote I considerando os erros independentes ( o que não ocorre!!) mod1 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas, start=c(a=1.1179, b=0.0562, c=-9.5923, d=0.367, e=0.0226)) summary(mod1) # quadro de estimativas # função para o cálculo do R2 R2 <- function(residuals, observed){ 1-(sum(residuals)^2)/sum((observed-mean(observed))^2) } R2(residuals(mod1), gas$y) residuosnI<-summary(mod1)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod1) hist(residuosnI) plot( ACF(mod1, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - AR() mod2 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corAR1(), start=c(a=0.78, b=0.06, c=-5.5923, d=0.367, e=0.0226)) summary(mod2) # quadro de estimativas R2(residuals(mod2), gas$y) residuosnI<-summary(mod2)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod2) hist(residuosnI) plot(x,mod2$residuals) plot(mod2) plot( ACF(mod2, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - AR() mod3 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corAR1(form=~x), start=c(a=0.78, b=0.06, c=-5.5923, d=0.367, e=0.0226)) summary(mod3) # quadro de estimativas R2(residuals(mod3), gas$y) residuosnI<-summary(mod3)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod3) hist(residuosnI) plot(x,mod3$residuals) plot( ACF(mod3, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - AR2() mod4 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corARMA(p=2,q=0), start=c(a=0.78, b=0.06, c=-5.5923, d=0.367, e=0.0226)) summary(mod4) # quadro de estimativas R2(residuals(mod4), gas$y) residuosnI<-summary(mod4)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod4) hist(residuosnI) plot( ACF(mod4, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - ARIMA(p=2,q=1) mod5 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corARMA(p=2,q=1), start=c(a=1.1179, b=0.0562, c=-9.5923, d=0.367, e=0.0226)) summary(mod5) # quadro de estimativas R2(residuals(mod5), gas$y) residuosnI<-summary(mod5)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod5) hist(residuosnI) plot( ACF(mod5, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - ARIMA(p=2,q=2) mod6 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corARMA(p=2,q=2), start=c(a=1.1179, b=0.0562, c=-9.5923, d=0.367, e=0.0226)) summary(mod6) # quadro de estimativas R2(residuals(mod6), gas$y) residuosnI<-summary(mod6)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod6) hist(residuosnI) plot( ACF(mod6, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - ARIMA(p=3,q=1) mod7 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corARMA(p=3,q=1), start=c(a=1.1179, b=0.0562, c=-9.5923, d=0.367, e=0.0226)) summary(mod7) # quadro de estimativas R2(residuals(mod7), gas$y) residuosnI<-summary(mod7)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod6) hist(residuosnI) plot( ACF(mod7, maxLag =24), alpha = 0.05 ) # ------------------------------------------------------------------------ # Ajuste do modelo Logístico 2 pool considerando - AR() a ordem das obsevações mod8 <- gnls(y ~ a/(1+exp(2+4*b*(c-x)))+d/(1+exp(2+4*e*(c-x))), data=gas,correlation=corAR1(form=~1), start=c(a=0.78, b=0.06, c=-5.5923, d=0.367, e=0.0226)) summary(mod8) # quadro de estimativas R2(residuals(mod8), gas$y) residuosnI<-summary(mod8)$resid shapiro.test(residuosnI) ks.test(rnorm(577),residuosnI) plot(mod6) hist(residuosnI) -- Adriele Giaretta Biase. Mestre em Estatística e Experimentação Agropecuária - UFLA. Doutoranda em Estatística e Experimentação Agronômica - ESALQ/ USP Contato: (19) 8861-0619.
participantes (1)
-
Adriele Giaretta Biase