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

[Fernando Antonio de souza]

Estou no celular e não tenho como avaliar de forma mais completa. Mas
realizei  análise dados semelhante s a um tempo atrás e a estrutura
corAR1(form=~1|grupo) me atendeu, onde grupo é o fator de agrupamento (
acredito ser amostra...) uma vez que as medidas tomadas em tempos
diferentes tem em comum a mesma amostra.

Assim q puder detalho melhor o q fiz, se alguém não esclarecer antes
Em 25/10/2015 9:32 PM, "Adriele Giaretta Biase" <adrielegbiase em gmail.com>
escreveu:

> Olá,
>
>
> eu gostaria de saber se é possível usar estruturas de correlação nos
> resíduos quando ajusto modelos não-lineares em mensurações de dados de
> produção de gases *in vitro* (caracterizado por mensurações acumuladas ao
> longo do tempo) e sem repetição.
>
>
> Breve descrição do banco de dados de gases *in vitro*: o experimento foi
> conduzido em frascos de (125 mL)  que contem fruído ruminal com feno de
> alfalfa. Esses frascos são acoplados aos sensores que medem pressão dos
> gases no interior dos frascos com o conversor digital que registra
> automaticamente a cada 5 minutos, durante o período de 48 horas. As
> mensurações dos gases são, portanto acumuladas ao longo do tempo. Tenho
> quatro compostos (amostras) diferentes que chamei de (y, y2, y3 e y4), em
> que cada uma amostra contem 577 mensurações ao longo do tempo (x) que é
> dado em horas.
>
>
> Eu preciso ajustar modelos não-lineares  para cada amostra (exponencial,
> logístico e logístico bicompartimental). O modelo logístico
> bicompartimental foi o que melhor se ajustou sem usar estrutruras de
> correlação, porém, a análise residual possuí uma forte dependência (visto
> que os dados são acumulados, ou seja, a mensuração no tempo (t+1) depende
> da mensuração anterior do tempo t.
>
> Então, eu teria que entrar com estruturas de correlações no modelo
> logístico bicompartimental para melhorar a análise residual?
>
>
>  Eu tentei algumas estruturas usando o argumento “correlation”, da função
> gnls. No quadro abaixo segue o seguinte resumo das estruturas inseridas no
> modelo, com os respectivos resultados gerados pela FAC ou os erros
> retornados pelo software R.
>
> modelo
>
> correlation
>
> Erros
>
> FAC - residual
>
> M1
>
> sem estrutrura de correlação
>
> Sem erro
>
> Tendência : possivelmente independência e heterogeneidade
>
> M2
>
> AR1: corAR1()
>
> Erro: step halving factor reduced below minimum in NLS step
>
> -
>
> M3
>
> AR1: corAR1(form=~x)
>
> Sem erro
>
> Tendência : possivelmente independência e heterogeneidade
>
> M4
>
> AR2: corAR1(p=2, q=0)
>
> Coefficient matrix not invertible
>
> -
>
> M5
>
> ARIMA: corARMA(p=2,q=1)
>
> step halving factor reduced below minimum in NLS step
>
> -
>
> M6
>
> ARIMA: corARMA(p=2,q=2)
>
> Coefficient matrix not invertible
>
> -
>
> M7
>
> ARIMA: corARMA(p=3,q=1)
>
> Coefficient matrix not invertible
>
> -
>
> M8
>
> AR(): corAR1(form=~1)
>
> Erro: step halving factor reduced below minimum in NLS step
>
> -
>
>
>
> Observação:  Eu fiz a análise usando a metodologia de séries temporais
> (modelo linear não estacionário – com uma diferença d=1) e a FAC e FACP
> mostraram que os resíduos eram ruído branco, como nos  modelos ARIMA(p=1,
> d=1, q=0),  ARIMA(p=2, d=1, q=0) e ARIMA(p=3, d=1, q=0). Seria possível
> ajustar modelos não-lineares (logístico bicompartimental) usando estruturas
> de correlações como o ARIMA?
>
> Os dados e o script seguem em anexo com todos os modelos que tentei
> ajustar (somente para uma amostra, ou seja, correspondente as colunas
> chamadas de y e x). Agradeço se puderem me dar uma orientação.
>
>
>
>
>
> rm(list=ls(all=TRUE))
>
> #pacotes necessários para a análise.
>
> require(manipulate)
>
> require(nlstools)
>
> require(car)
>
> require(nlme)
>
> require(lmtest)
>
> # modelos estudados
>
> #a * (1 - exp(-b * (t - c))) # modelo exponencial
>
> #a /(1+exp(2+4*b*(c-t)))     # modelo logistico
>
> #a/(1+exp(2+4*b*(c-t)))+d/(1+exp(2+4*e*(c-t)))  # modelo logistico 2 pool
>
>
>
> #-----------------------------------------------------------------------
>
> setwd('E:\\autoregressivo_modelo')
>
> gas01<-read.csv("gas.csv",sep = ",", header = TRUE)
>
> gas01
>
> head(gas01)
>
> attach(gas01)
>
> t=Time/60    # para deixar o tempo em horas
>
> dim(gas01)
>
>
>
> gas=cbind(gas01,t)
>
> names(gas) <- c("Time","y","y2","y3", "y4", "x") # variáveis que quero
> ajustar y (variavel dependente gases) e x o tempo em horas
>
> head(gas)
>
> dim(gas)
>
> attach(gas)
>
>
>
> #===================================================================
>
>
>
> # Ajuste usando o gnls
>
> #Ajuste do modelo exponencial -lote I considerando os erros independentes
>
> 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
>
> R2(residuals(mod1), gas$y)
>
> residuosnI<-summary(mod1)$resid
>
> shapiro.test(residuosnI)
>
> ks.test(rnorm(577),residuosnI)
>
> plot(mod1)
>
> hist(residuosnI)
>
>
>
>
>
> # 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)
>
> # ------------------------------------------------------------------------
>
>
>
>
>
> # 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 )
>
>
>
>
>
> Agradeço desde já,
>
> Adriele.
>
> --
> 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.
>
> _______________________________________________
> R-br mailing list
> R-br em listas.c3sl.ufpr.br
> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
> código mínimo reproduzível.
>
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20151025/7d8df2f8/attachment-0001.html>