[R-br] tratamento de pontos discrepantes (pesos de bovinos)
Adriele Giaretta Biase
adrielegbiase em gmail.com
Quinta Julho 20 16:42:59 -03 2017
Boa tarde pessoal,
estou trabalhando com um banco de dados que possuem várias coletas de pesos
ao longo do tempo de vários bois.
O meu objetivo é deixar um ajuste simples (pode até ser um modelo não-
linear ao invés de uma regressão simples), porém, preciso remover os
possíveis pontos discrepantes que atrapalham o ajuste, de forma recursiva
(automática). O algoritmo será executado em vários momentos, enquanto
acompanha o crescimento do animal. O mais importante seria a remoção dos
*outiliers* do banco de dados. Eu criei um algoritmo para fazer isso.
Gostaria de saber se alguém possui uma ideia melhor ou sugestão. Segue a
função criada.
*OBS:* Muitos pontos discrepantes ainda ficam no banco, isso devido a um
erro operacional da balança que pesa os animais. Precisaria de algo mais
robusto para remover os pontos destoantes de forma recursiva para todos os
animais.
DADOS é um data.frame com as seguintes colunas: BRINCO (identificação ou
numero do animal), DATA (dia da pesagem do peso do animal) e PESO (peso do
animal correspondente ao dia)
fit_linear <- function(DADOS){
fit_aux = c()
BRINCO = as.numeric(as.vector(levels( as.factor(DADOS[,"BRINCO"]))))
pdf(file = "ajuste%03d.pdf") # salvar os gráficos em pdf
for (i in BRINCO){
dados = subset(DADOS, BRINCO == i, selec = c( BRINCO, DATA, PESO))
j=1
while (j<=30){
dias = as.vector(as.numeric(as.factor(dados[,"DATA"])))
pesos = as.numeric(as.vector(dados$PESO))
dad = cbind(dias, pesos)
mod <- lm(pesos~dias)
INFLUENTES = influence.measures(mod)
dados_analise = cbind(INFLUENTES$is.inf, dad)
dados_influentes = as.data.frame(dados_analise)
colnames(dados_influentes) = c("x1", "x2", "x3", "x4", "x5", "x6", "DATA",
"PESO")
dados = subset(dados_influentes, (x1==FALSE & x2==FALSE & x3==FALSE &
x4==FALSE & x5==FALSE & x6==FALSE),
selec = c(x1, x2, x3, x4, x5, x6, DATA, PESO))
j = j+1
}
res <- residuals(mod)
h <- hatvalues(mod)
P = length(mod$coefficient); N=length(peso); P; N; hc<-3*P/N;hc
limt <- list(c(DFB=2/sqrt(N),FDFits=2*sqrt(P/N),cov.r=3*P/N,
Cook=qf(0.5,2,8, lower.tail = TRUE, log.p = FALSE),hat=3*P/N))
N=length(peso); P=length(mod$coefficient)
rs <- rstudent(mod)
h <- lm.influence(mod)$hat; lc <- 3*P/N
minrs=min(min(rs),-3)
maxrs=max(max(rs),3)
ymin=minrs-.1
ymax=maxrs+.1
maxh=max(max(h),lc)
minh=min(h)
xmin=minh-.1
xmax=maxh+.1
par(mfrow=c(1,1))
plot(c(xmin,xmax),c(ymin,ymax), type="n", xlab="h - leverage",
ylab="RStudent", main= i)
abline(h=-2.5, col="red")
abline(h=2.5,col="red"); abline(v=lc, col="blue")
points(h,rs)
plot(dad, main= i)
abline(mod,lty=2)
}
dev.off()
}
Agradeço,
--
Adriele Giaretta Biase.
Mestre em Estatística e Experimentação Agropecuária - UFLA.
Doutora em Estatística e Experimentação Agronômica - ESALQ/ USP
Contato: (19) 98861-0619.
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20170720/d9a94728/attachment.html>
Mais detalhes sobre a lista de discussão R-br